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Abstract 


A  thermal  conductive  heating  (TCH)  pilot  test  was  conducted  at  the  Naval  Air  Warfare  Center 
(NAWC)  in  West  Trenton,  New  Jersey  in  2009  in  collaboration  with  TerraTherm,  Inc.,  the  Naval 
Facilities  Engineering  Services  Center  and  the  United  States  Geological  Survey.  The  NAWC  site 
was  historically  used  as  a  jet  engine  testing  facility  from  the  mid-1950s  to  the  late  1990s.  During 
this  time,  the  subsurface  was  contaminated  with  trichloroethylene  (TCE)  which  was  a  common 
solvent  used  at  the  facility.  The  pilot  test  consisted  of  15  heater/extraction  wells  installed  to  a 
depth  of  16.8  m  in  weathered  mudstone  and  operated  for  102  days.  Rock  core  samples  were  taken 
pre-  and  post-remediation  to  measure  the  initial  TCE  concentrations  and  evaluate  the  effect  the 
TCH  pilot  test  had.  The  data  collected  during  the  pilot  test  was  used  to  create  a  two-dimensional 
(2D)  finite  difference  model  using  TMVOC.  TMVOC  is  part  of  the  TOUGH  2  family  of  codes 
and  is  a  numerical  model  that  is  capable  of  simulating  multiphase  flow,  heat  transfer  and  transport 
of  volatile  organic  compounds  in  three-dimensional  heterogenous  porous  media  or  fractured  rock. 

The  2D  model  was  used  as  a  screening  model  to  investigate  TCE  removal  from  the  rock  matrix 
when  heating  for  100  days  with  a  similar  heating  pattern  to  what  was  employed  at  the  NAWC 
site.  The  numerical  domain  incorporated  three  primary  fractures  with  competent  bedrock  in 
between.  As  the  test  pilot  was  conducted  in  the  weathered  bedrock  zone,  a  sensitivity  analysis 
was  first  completed  on  the  matrix  permeability  to  help  to  match  the  TCE  removal  from  the  pilot 
test.  The  pilot  test  had  a  63.5%  removal  of  TCE  from  the  study  area  compared  to  67%  from  the 
baseline  model.  A  limited  sensitivity  analysis  was  completed  which  investigated  how  the  matrix 
porosity  and  rate  of  energy  application  would  have  on  the  success  of  TCE  removal  from  the  rock 
matrix.  It  revealed  that  the  TCE  removal  increases  with  increased  matrix  porosity  and  increased 
rate  of  energy  application. 
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Chapter  1 


Introduction 

Groundwater  and  soil  contamination  by  dense,  non-aqueous  phase  liquids  (DNAPLs)  has  been  a 
continuing  concern  in  North  America  since  the  early  1980s.  Historical  disposal  methods, 
regulations  and  governance  have  resulted  in  tens  of  thousands  of  sites  impacted  by 
polychlorinated  biphenyls  (PCBs),  creosote,  coal  tar,  and  various  chlorinated  solvents.  The  focus 
of  this  research  is  on  the  remediation  of  the  chlorinated  solvent  trichloroethylene  (TCE)  in 
fractured  bedrock.  A  numerical  model  capable  of  simulating  three-phase  non-isothermal  flow  of 
water,  gas  and  multi-component  mixtures  of  volatile  organic  compounds  (VOCs)  in 
multidimensional  heterogeneous  porous  media  was  used  to  model  a  Thermal  Conductive  Heating 
(TCH)  field  pilot  test  that  was  conducted  at  a  contaminated  fractured  bedrock  site.  The  pilot  test 
took  place  at  the  Naval  Air  Warfare  Centre  (NAWC)  in  West  Trenton,  NJ  which  was  previously 
used  as  a  jet  engine  testing  facility,  and  decommissioned  in  1998.  The  results  from  the  model 
have  been  compared  to  the  field  test  results  to  ascertain  what  properties  of  the  bedrock  and  what 
features  of  the  technology  are  most  influential  with  respect  to  mass  and  concentration  reduction. 

Organic  liquids  immiscible  with  water  are  referred  to  as  non-aqueous  phase  liquids  (NAPLs).  If 
the  NAPL  is  less  dense  than  water,  it  is  referred  to  as  a  light,  non-aqueous  phase  liquid  (LNAPL); 
if  it  is  denser  than  water,  it  is  referred  to  as  a  dense,  non-aqueous  phase  liquid  (DNAPL).  There 
are  several  classes  of  DNAPLs  which  include  chlorinated  solvents,  coal  tar,  creosote,  and  PCBs 
(Kueper  et  al.,  2004).  In  addition  to  being  denser  than  groundwater,  chlorinated  solvent  DNAPLs 
often  have  a  lower  viscosity  than  water,  resulting  in  potentially  rapid  migration  rates  below  the 
water  table,  particularly  in  fractured  rock  (Kueper  and  McWorter,  1991).  It  is  difficult  to  predict 
the  migration  pathways  of  DNAPL  when  the  geology  at  a  site  is  complex  as  it  often  is  with 
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fractured  rock  (Domenico  &  Schwartz,  1 990).  Upon  release  to  the  subsurface,  DNAPL  will  come 
to  rest  as  both  disconnected  blobs  and  ganglia  of  organic  liquid  referred  to  as  residual  DNAPL, 
and  in  continuous  disfributions  referred  to  as  pooled  DNAPL  (Kueper  et  al,  1997).  Both  residual 
and  pooled  DNAPL  slowly  dissolve  into  flowing  groundwater,  giving  rise  to  aqueous  phase 
plumes. 

Due  to  inadequate  knowledge  and  technology  in  the  20th  century,  it  was  unknown  that  chlorinated 
solvents  could  result  in  long-term  contamination  of  the  subsurface.  Poor  disposal  and  storage 
procedures  resulted  in  many  releases  to  the  subsurface  that  were  not  properly  remediated  (Kueper 
et  al,  2004).  There  are  many  sites  in  the  United  States  and  Canada  that  are  contaminated  with 
chlorinated  solvent  DNAPLs  and  are  in  need  of  remediation. 

TCE  is  a  common  industrial  degreaser  and  dry  cleaning  fluid  that  has  been  widely  used  by  a 
variety  of  industries  since  the  1940s  (Kueper  et  al,  2004).  TCE  has  been  shown  to  be  harmful  to 
humans  if  repeated  or  high  levels  of  exposure  occur  (U.S.  Department  of  Health  and  Human 
Services,  1997).  Exposure  can  lead  to  acute  and  chronic  health  effects  (U.S.  Department  of 
Health  and  Human  Services,  1997).  Stringent  drinking  water  standards  have  been  put  in  place  for 
allowable  levels  of  TCE,  due  to  potential  short-term  and  long-term  health  effects.  In  the  province 
of  Ontario,  for  example,  the  maximum  allowable  concentration  of  TCE  in  potable  water  is  5  pg/L 
(MOE,  2006,  Table  2). 

TCE  was  used  at  the  NAWC  facility  from  the  mid-1950s  to  the  late  1990s  where  jet  engine 
testing  occurred  (Lacombe,  2000).  The  activities  at  the  site  resulted  in  TCE  and  other  chemicals 
being  released  into  the  subsurface.  A  pump-and-treat  system  has  been  in  place  at  the  NAWC  site 
since  mid- 1990  (Lacombe,  2011).  The  site  is  underlain  by  fractured  sedimentary  bedrock, 


primarily  mudstone.  Mudstone  is  a  sedimentary  rock  that  is  composed  of  lithified  mud  and  shale 
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(Aplin  et  al.,  1999).  A  fissile  mudstone  is  known  as  a  mud-shale,  while  a  non-fissile  mudstone  is 
referred  to  as  mudstone  (Aplin  et  al.,  1999).  Mudstones  are  classified  primarily  on  grain-size  and 
are  typically  composed  of  1/3  to  2/3  silt  sized  fractions. 

There  are  several  remediation  methods  that  have  been  used  to  address  DNAPLs  in  the  subsurface. 
These  include  pump-and-treat,  hydraulic  displacement,  surfactant  and  co-solvent  flushing, 
enhanced  in  situ  bioremediation  (EISB),  in-situ  chemical  oxidation  (1SCO),  steam-enhanced 
extraction  (SEE),  electrical  resistance  heating  (ERH),  and  thermal  conductive  heating  (TCH).  The 
remedial  success  of  many  of  these  methods  in  bedrock  relies  on  the  ability  to  flush  fluids  through 
the  fracture  network  (Pankow  and  Cherry,  1 996). 

SEE,  ERH  and  TCH  are  remediation  methods  that  use  heat  to  increase  the  subsurface 
temperature,  thereby  enhancing  DNAPL  removal  (Triplett  Kingston  et  al,  2010).  In  the 
application  of  SEE,  a  network  of  injection  and  extraction  wells  are  used  whereby  steam  is 
injected  into  the  subsurface  and  the  vapors  and  liquids  are  recovered  for  treatment  (Johnson  et  al, 
2009).  The  process  creates  a  steam  zone  in  the  treatment  area  (Nilsson  et  al,  2011).  The 
technology  relies  on  the  ability  to  deliver  steam  to  the  contaminated  areas  and  therefore  may  be 
difficult  to  apply  in  fractured  rock.  A  general  limiting  factor  of  SEE  is  that  the  boiling  point  of 
water  is  the  maximum  temperature  that  can  be  obtained. 

ERH  is  similar  to  SEE  in  that  the  maximum  obtainable  temperature  in  the  subsurface  is  the 

boiling  point  of  water.  This  technology  uses  a  network  of  electrodes  to  conduct  electrical  current 

through  the  subsurface,  along  with  a  network  of  extraction  wells  to  collect  liquids  and  vapours 

(Johnson  et  al,  2009).  This  technology  performs  best  when  there  is  ample  moisture  in  the 

subsurface  to  maintain  electrical  conductivity.  The  effectiveness  of  ERH  can  be  adversely 

affected  by  groundwater  flow.  If  the  groundwater  velocity  is  sufficiently  high,  heat  will  be 
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carried  away  from  the  treatment  zone  at  a  rate  sufficient  to  inhibit  boiling.  ERH  is  generally  not 
applicable  to  fractured  rock  because  of  the  excessively  high  electrical  resistivity  of  most  rock 
types. 

TCH  utilizes  a  network  of  wells  containing  a  heating  coil  to  heat  the  subsurface  through 
conduction  to  temperatures  reaching  500  degrees  Celsius  (°C)  (Heron  et  al,  2009).  This  process  is 
highly  effective  for  volatile  organic  compounds  such  as  chlorinated  solvents  as  the  dominant 
removal  mechanism  is  by  vaporization  of  the  chlorinated  solvent  and  collection  of  the  subsequent 
vapors  by  extraction  wells.  Similar  to  ERH,  the  success  of  thermal  conductive  heating  is  reduced 
with  increased  groundwater  flow.  TCH  has  been  successfully  implemented  at  several 
contaminated  sites.  For  example,  in  Memphis,  Tennessee  a  site  had  levels  of  chlorinated 
compounds  at  concentrations  greater  than  1000  mg/kg.  The  site  was  underlain  by  low- 
permeability  loess  which  was  in  turn  underlain  by  high-permeability  alluvium  composed  of 
sands,  sands  and  gravels,  and  fluvial  deposits  (Heron  et  al,  2009).  A  volume  of  38,200  m3  was 
treated  using  TCH  for  177  days  (Heron  et  al,  2009).  The  site  was  remediated  to  concentrations 
less  than  0.01  mg/kg,  representing  a  99.99%  reduction  in  concentration. 

1.1  Research  Objectives 

Thermal  conductive  heating  applications  have  been  implemented  on  many  occasions  in  porous 
media  with  successful  results  (Johnson  et  al.,  2009),  but  there  is  relatively  little  experience  in 
fractured  rock.  .  To  evaluate  the  performance  of  TCH  in  fractured  rock,  a  field  pilot  test  was 
earned  out  by  TerraTherm,  Inc.  in  2009  in  collaboration  with  the  Naval  Facilities  Engineering 
Services  Center,  Queen’s  University  and  the  United  States  Geological  Survey  at  the  former  Naval 
Air  Warfare  Center  (NAWC)  in  West  Trenton,  New  Jersey.  Queen’s  University  (David 
Rodriguez  and  Ashley  McKenzie)  assisted  in  data  collection  by  measuring  the  trichloroethylene 
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concentrations  in  the  rock  pre-  and  post-remediation  using  methods  developed  by  Rodriguez  and 
Kueper,  2012).  Using  data  gathered  from  site  characterization  activities  at  the  NAWC  site  and 
the  TCH  pilot  test,  a  finite  difference  model  was  created  using  TMVOC  to  evaluate  the  reduction 
of  TCE  concentrations  in  the  rock  matrix  upon  application  of  TCH. 

TMVOC  is  a  numerical  model  that  is  capable  of  simulating  the  flow  and  transport  of  multi- 
component  VOCs  in  all  three  phases  (gas-water-NAPL),  above  and  below  the  water  table  in 
three-dimensional  heterogeneous  porous  media,  or  fractured  bedrock.  It  is  able  to  simulate  non¬ 
condensable  gases  as  well  as  advection,  dispersion,  soiption,  phase  partitioning  and 
biodegradation  (Pruess  and  Battistelli,  2002).  TMVOC  allows  for  irregular  grids  to  be  created  by 
using  the  Integral  Finite  Difference  Method  for  space  discretization.  The  model  can  simulate  heat 
transfer  applications  such  as  thermal  conductive  heating  and  steam-enhanced  extraction 
remediation  technologies. 

A  radial  symmetric  screening  model  of  the  pilot  test  completed  at  the  NAWC  site  was  created. 
The  temperature  and  mass  removal  data  collected  during  the  test  pilot  were  used  as  input 
parameters  to  create  the  baseline  model  conditions.  A  limited  sensitivity  analysis  was  completed 
on  the  2D  model  to  evaluate  which  parameters  are  most  influential  on  TCE  removal  from  the  rock 
matrix.  The  varied  parameters  in  the  sensitivity  analysis  were  matrix  permeability,  matrix 
porosity  and  energy  rate.  From  the  results  of  this  study,  recommendations  for  future  TCH 
remediation  in  fractured  bedrock  aquifers  have  been  determined  as  well  as  future  modeling 
recommendations. 
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1.2  Organization 


Organization  of  the  remaining  Chapters  and  Appendices  are  as  follows: 

•  Chapter  2  is  a  literature  review  that  provides  background  information  pertaining  to  the 
work  completed  for  this  project. 

•  Chapter  3  has  been  written  as  a  standalone  manuscript  form  intended  for  publication  in 
Ground  Water  Monitoring  and  Remediation.  Ashley  McKenzie  is  the  primary  author  of 
this  paper  and  Bernard  Kueper  and  Ronald  Falta  are  the  co-authors. 

Supporting  information  for  the  project  as  well  as  the  manuscript  is  contained  in  the  Appendices 
which  include: 

•  Appendix  A  -  TMVOC  Governing  Equations  which  provides  a  brief  summary  of  the 
equations  that  are  used  by  TMVOC; 

•  Appendix  B  -  Model  Testing  which  examined  if  the  numerical  model  was  calibrated 
correctly  by  recreating  an  existing  TMVOC  model  from  the  Users  Manual; 

•  Appendix  C  -  TMVOC  Model  Input  Property  Definitions  which  further  explain  the  input 
parameters  required  for  the  numerical  model  simulations; 

•  Appendix  D  -  Calculations  that  were  required  for  select  model  input  parameters;  and, 

•  Appendix  E  -  Fixed  versus  No-Flow  Bottom  Boundary  Condition  which  examines  the 
two  different  bottom  boundary  conditions  considered  when  calibrating  the  baseline 
simulation. 
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Chapter  2 


Literature  Review 

Groundwater  and  soil  contamination  is  a  growing  issue  in  North  America.  Previous  disposal 
methods,  regulations  and  governance  have  led  to  many  types  of  contamination  in  the  subsurface. 
The  focus  of  this  research  is  on  chlorinated  solvents,  which  have  been  used  as  industrial 
degreasers  since  the  early  1900s.  Remediation  of  these  solvents  can  be  challenging  and  is  very 
site  specific.  Thermal  conductive  heating  (TCH)  is  a  technology  that  has  successfully  removed 
chlorinated  solvents  from  porous  media  sites.  It  is  a  method  by  which  a  grid  of  wells  containing  a 
heating  element  heats  the  subsurface  to  above  the  boiling  temperature  of  the  chemical,  with  an 
extraction  system  at  the  surface  to  collect  the  vapors  for  further  treatment.  It  is  important  to 
understand  the  characteristics  of  chlorinated  solvents,  thermal  and  fractured  rock  properties,  as 
well  as  TCH  before  implementing  this  technology. 

2.1  Chlorinated  Solvents 

Chlorinated  solvents  were  introduced  to  North  America  in  the  20th  century  when  they  were  used 
as  degreasers  in  manufacturing  plants,  and  for  many  products  such  as  fire  extinguishers  and 
refrigeration  uses  (Jackson  &  Dwarakanath,  1999).  Chlorinated  solvents  are  commonly  referred 
to  as  dense  non  aqueous  phase  liquids  (DNAPLs)  that  are  composed  of  chlorinated  hydrocarbons 
and  can  be  single  or  multi-component  compounds  (Kueper  et  al,  2004).  The  most  common 
chlorinated  solvents  are  trichloroethylene  (TCE),  perchloroethylene  (PCE)  and  carbon 
tetrachloride  (CT),  which  were  used  since  the  1960s  by  industry  for  a  variety  of  puiposes 
including  metal  degreasing  and  dry  cleaning.  Due  to  inadequate  knowledge  and  technology  in 
the  20th  century,  it  was  unknown  that  chlorinated  solvents  could  contaminate  the  subsurface, 
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therefore  disposal  and  storage  procedures  were  poor,  and  many  accidental  spills  occurred  which 
were  not  properly  cleaned  up  (Kueper  et  al,  2004). 

Chlorinated  solvents  are  heavier  than  water,  which  has  a  density  of  998  kg/m3  and  a  viscosity  of 
0.01005  g/cm-s  at  20  degrees  Celsius,  as  compared  to  common  chlorinated  solvent  DNAPLs  that 
have  densities  ranging  from  1100  kg/nr3  to  1600  kg/m3,  and  viscosities  ranging  from  0.0057 
g/cm-s  to  0.01  g/cm-s  respectively  (Fetter,  2001;  Kueper  et  al,  2004).  With  a  higher  density  than 
water,  chlorinated  solvents  are  able  to  migrate  below  the  watertable,  with  potentially  rapid  rates 
of  migration  in  the  subsurface  (Kueper  et  al,  2004). 


There  are  many  sites  in  the  United  States  and  Canada  that  are  contaminated  with  DNAPLs,  and 
are  in  need  of  remediation.  It  is  difficult  to  predict  the  migration  pathways  of  a  DNAPL  when  the 
geology  at  the  site  is  complex  (Domenico  &  Schwartz,  1990).  Figure  2-1  is  an  example  of  how 
DNAPL  migration  may  occur  in  the  subsurface  (Kueper  and  McWhorter,  1991). 


ON  API.  RELEASE 


Figure  2-1  -  Migration  pathways  due  to  DNAPL  release  into  the  subsurface  (Kueper,  & 

McWhorter,  1991). 
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Residual  and  pooled  DNAPL  will  form  as  the  DNAPL  migrates  through  the  subsurface  (Kueper 
et  al.,  1997).  Residual  DNAPL  is  formed  when  the  DNAPL  migrates  into  pore  spaces  or  fractures 
and  becomes  trapped  by  capillary  forces  as  isolated  blobs  and  ganglia  of  liquid.  Residual  DNAPL 
can  occur  in  the  saturated  or  unsaturated  zone,  and  the  amount  of  residual  created  varies 
according  to  many  factors  such  as  interfacial  tension,  the  viscosity  and  density  of  the  DNAPL, 
and  the  permeability  of  the  medium  (Kueper  et  al.,  1997).  Pooled  DNAPL  represents  a 
continuous  distribution  of  liquid  through  the  pore  spaces,  and  is  found  above  capillary  barriers 
typically  associated  with  variations  in  permeability  (Kueper  et  al.,  1997). 

Potable  groundwater  for  generic  site  conditions  have  maximum  acceptable  concentrations  of  0.5 
|Jg/L  for  TCE  and  PCE,  and  0.2  (Jg/L  for  CT  as  set  by  the  Ministry  of  Environment  (MOE,  2011, 
Table  1).  Repeated  exposure  to  these  chemicals  above  the  acceptable  concentrations  can  lead  to 
acute  and  chronic  health  effects  which  may  even  lead  to  death  (Ajo-Franklin  et  al.,  2006). 

The  focus  of  this  study  is  on  thermal  remediation  of  TCE.  The  production  of  TCE  began  in  1923 
for  use  in  adhesives  and  for  metal  cleaning  and  degreasing  (Pankow  &  Cherry,  1996).  TCE  is  a 
non-flammable,  colorless  liquid  with  a  sweet  odor  and  burning  taste  (Irwin,  1997).  In  the  1970s 
production  of  TCE  declined  when  several  health  concerns  arose  (Pankow  &  Cherry,  1996).  TCE 
can  break  down  into  other  dangerous  substances  such  as  dichloroethylene,  chloromethane,  and 
vinyl  chloride,  which  are  also  dangerous  if  exposure  occurs,  and  causes  further  soil  contamination 
(Irwin,  1997).  Table  2.1  displays  the  properties  of  TCE  compared  to  water.  When  remediating 
TCE,  these  chemical  properties  are  important  to  consider  as  they  may  be  affected  by  temperature 
and  pressure  changes. 
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Table  2-1  -  Physical  and  chemical  properties  of  trichloroethylene  in  comparison  to  water. 


Parameter 

TCE 

Water 

Chemical  Formula 

C2HC13 

H20 

Liquid  density  at  20°C 

1460  kg/m3  (a) 

998  kg/m3  (a'd) 

Viscosity 

0.0057  g/cm-s(b’c) 

0.010050  g/crn-s  <d) 

Aqueous  solubility 

1100  mg/L  (c) 

Not  defined 

Surface  tension  with  air 

29.5  dynes/cm <a) 

72.0  dynes/cm (a) 

Maximum  contaminant  level 

0.005  mg/L  (a) 

Not  defined 

Ajo-Franklin,  Geller,  &  Harris  (2006). 
b  Jackson,  &  Dwarakanath  (1999). 

L  Kueper  et  al  (2004). 
d  Fetter  (2001). 

2.2  Fluid  Properties 

2.2.1  Viscosity 

Viscosity  is  the  measure  of  the  resistance  of  a  fluid  to  flow,  which  is  caused  by  internal  friction 
within  a  fluid  (Knox  et  al.,  1993,  Mercer  and  Cohen,  1990).  Chlorinated  solvents  are  less  viscous 
than  water  with  a  viscosity  ranging  from  0.57  to  1.0  cP  compared  to  that  of  water  at  1.0  cP 
(Kueper  et  al.,  2004).  A  lower  viscosity  results  in  the  chemical  penetrating  the  subsurface  more 
rapidly  with  easier  migration  (Mercer  and  Cohen,  1990;  Pankow  and  Cherry,  1996).  The  viscosity 
of  a  fluid  is  inversely  proportional  to  temperature.  As  the  temperature  of  DNAPL  increases,  the 
viscosity  decreases  which  results  in  a  less  viscous  and  more  mobile  fluid  (O’Carroll  &  Sleep, 
2007). 

2.2.2  Interfacial  Tension  and  Wettability 

Interfacial  tension  and  wettability  are  two  properties  of  chlorinated  solvents  that  are  closely 
related.  Interfacial  tension  is  the  force  formed  at  the  interface  between  two  immiscible  liquids, 
which  prevents  complete  mixture  of  the  fluids  (Mercer  and  Cohen,  1990).  Typical  values  for 
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interfacial  tension  for  DNAPLs  and  water  existing  in  the  subsurface  are  between  5  and  35  mN/m 
(Kueper  et  al.,  1997).  These  are  less  than  textbook  values  because  of  the  fact  that  small  amounts 
of  impurities  typically  lead  to  a  reduction  in  interfacial  tension  (Dou  et  al.,  2008).  The  interfacial 
tension  between  air  and  water  is  72  mN/m,  which  is  significantly  greater  than  that  between 
DNAPLs  and  water.  The  interfacial  tension  between  two  fluids  will  also  be  affected  by  the 
temperature.  For  example,  as  the  temperature  increases,  the  surface  tension  of  air-water  decreases 
(Vargaftik  et  al.,  1983).  The  wettability  is  dependent  on  the  interfacial  tension  of  a  chemical  and 
is  defined  as  the  ability  for  a  chemical  to  spread  on  a  surface  when  present  with  another  chemical 
(Dou  et  al,  2008;  Mercer  and  Cohen,  1990).  It  is  measured  by  the  contact  angle  between  the  two 
fluids  where  the  fluid-fluid  interface  meets  the  solid  surface  as  displayed  in  Figure  2-2. 


Figure  2-2  -  (top)  water-wet  contact  angle;  (bottom)  NAPL-wet  contact  angle  (modified 

from  Mercer  &  Cohen,  1990). 
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The  two  fluids  are  classified  as  either  wetting  or  non-wetting  depending  on  the  contact  angle. 
Young’s  equation  represents  the  relationship  between  the  contact  angle  and  interfacial  tension  as 
follows  (Mercer  and  Cohen,  1990): 


(CTns  °ws)  (2-1) 

COS0  =  - 

°nw 

where  (/>  is  the  contact  angle,  a  is  the  interfacial  tension,  n  represents  the  NAPL  phase,  s 
represents  the  solid  and  w  represents  the  water  phase.  When  the  contact  angle  is  less  than  70 
degrees,  the  fluid  through  which  that  angle  is  measured  is  considered  to  be  the  wetting  fluid,  and 
the  other  is  considered  to  be  the  non-wetting  fluid.  When  the  contact  angle  is  close  to  90  degrees, 
the  system  is  considered  to  be  neutral-wet  (Mercer  and  Cohen,  1990;  Donaldson  and  Alam, 
2008). 

2.2.3  Phase  Change 

When  a  liquid  is  subjected  to  a  change  in  temperature  or  pressure,  it  can  cause  the  liquid  to 
change  phases  into  a  vapor.  Temperature  and  pressure  are  very  important  parameters  when 
examining  DNAPLs  as  they  are  dependent  for  liquid-vapor  mixtures.  A  distinction  between  the 
liquid-vapor  phases  can  not  always  be  made  and  this  point  is  known  as  the  critical  point.  It  is  the 
exact  moment  when  the  change  occurs.  Figure  2-3  shows  three  different  phases:  solid,  liquid,  and 
vapor  for  a  pure  substance  when  it  is  subjected  to  pressure  and  temperature  in  relation  to  the 
volume. 
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Figure  2-3  -  Pressure-volume  diagram  for  a  pure  substance  (Smith  et  al,  1996). 


The  phase  changes  due  to  temperature  and  pressure  can  be  describe  by  the  equation  of  state  and  is 
described  by  Smith  et  al.  (1996)  as: 


dV  = 


(2.2) 


where  V  is  the  volume,  T  is  the  temperature,  and  P  is  the  pressure.  The  equation  of  state  can  be 
solved  for  pressure,  temperature,  or  volume  and  it  is  a  function  of  the  two  remaining  variables 
where  the  “partial  derivatives  have  definite  physical  meaning  and  are  measureable  quantities” 
(Smith  et  al.,  1996). 
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Boiling  and  condensation  is  an  important  consideration  in  phase  change.  The  heat  that  is  absorbed 
or  released  in  the  liquid-vapor  phase  is  described  by  Naterer  (2003)  as: 


(2.3) 


where  hfg  is  the  latent  heat  of  vaporization,  q  is  the  rate  of  heat  transfer  during  boiling,  and  m  is 
the  mass  of  liquid  boiled  per  unit  time  when  examining  the  vapor  phase  and  it  is  defined  as  the 
rate  of  vapor  condensation  on  a  cooled  surface  if  evaluating  the  condensation  phenomenon 
(Naterer,  2003). 

2.3  Heat  Transfer 
2.3.1  Conduction 

Conduction  is  a  form  of  heat  transfer  that  is  driven  by  the  temperature  gradient;  when  the  heat  (or 
energy)  moves  from  a  high  energy  region  to  a  low  energy  region  (Holman,  1997).  It  is  defined  by 
Fourier’s  Law  of  conduction  as  defined  by  Thomas  (1980)  as: 

q  =  -kAVT  (2.4) 


here  q  is  the  heat  transfer  rate  in  units  of  W/m-C,  k  is  the  thermal  conductivity  of  the  material,  A 
is  the  area  of  the  object  being  examined  and  VT  is  the  temperature  gradient.  The  second  law  of 
Thermodynamics  states  that  heat  is  transferred  in  the  direction  of  decreasing  temperature  or 
“downstream”  which  means  that  a  negative  sign  is  present  in  the  equation  to  make  the  heat 
transferrate  positive  (Thomas,  1980). 
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Equation  2.4  can  be  applied  to  a  three-dimensional  system  where  the  thermal  conductivity,  K,  is 
held  constant  to  become: 

q2t  q2T  d2T  q  _igT  (2.5) 

dx2  +  dy2  +  dz2  +  K  a  dr 

where  a  is  the  thermal  diffusivity  of  the  material.  An  increase  in  the  thermal  diffusivity  indicates 
that  the  temperature  will  diffuse  through  the  material  at  a  faster  rate  (Simons,  2007). 

2.3.2  Convection 

Convection  is  a  physical  process  of  heat  transfer  that  occurs  at  the  surface  of  a  moving  fluid  and 
can  only  occur  in  gases  or  liquids  (Thomas,  1980;  Simons,  2007).  The  equation  used  to  quantify 
convection  is  also  known  as  Newton’s  law  of  cooling  and  is  defined  by  Holman  (1997)  as: 


qc  =  hA(Tw  -  Tm)  (2.6) 

where  q  is  the  convective  heat  transfer  rate,  h  is  the  convection  heat  transfer  coefficient,  A  is  the 
surface  area,  Tw  is  the  temperature  of  the  wall  and  Tx  is  the  temperature  of  the  fluid.  Convection 
is  dependent  on  viscosity,  thermal  conductivity,  specific  heat  and  the  density  of  the  fluid 
(Holman,  1997).  There  are  two  types  of  convection:  forced  or  natural.  Forced  convection  is  when 
the  motion  or  flow  is  caused  by  mechanical  or  induced  external  forces,  and  natural  convection  is 
caused  by  buoyancy  forces  which  result  in  a  temperature  induced  density  gradient  within  the 
fluids  (Thomas,  1980;  Simons,  2007).  When  examining  thermal  remediation,  convection  would 
be  caused  by  the  migration  of  steam  produced  by  boiling  groundwater. 
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2.3.3  Radiation 


Objects  emit  and  absorb  electromagnetic  radiation,  which  occurs  due  to  the  vibrational  and 
rotational  molecular,  atomic  and  subatomic  particle  movements  (Thomas,  1980;  Simons,  2007). 
Thermal  radiation  is  emitted  by  physical  matter  (Thomas,  1980).  Radiation  can  only  occur 
through  transparent  or  partially  transparent  bodies  when  electromagnetic  radiation  is  involved 
(Simons,  2007). 

2.3.4  Thermal  Diffusivity 


Thermal  diffusivity  is  a  temperature  dependent  property  that  governs  the  transient  heat  transport 
equation  [2.5].  It  is  directly  related  to  the  thermal  conductivity  of  the  material,  k,  and  inversely 
related  to  the  specific  heat  at  a  constant  pressure,  cp.  It  is  expressed  by  Cermak  and  Rybach 
(1982)  as: 


cpP 


(2.7) 


where  a  is  the  thermal  diffusivity  and  p  is  the  density  of  the  material.  Studies  by  Durham  et  al. 
(1987)  have  shown  that  thermal  diffusivity  is  dependent  on  temperature  at  different  pressures.  It 
was  found  that  as  the  temperature  increases,  the  thermal  diffusivity  decreases.  The  experiment 
performed  by  Durham  et  al.  (1987)  also  found  that  the  diffusivity  increases  when  the  pressure  on 
the  rock  material  increases,  although  the  change  is  not  as  significant  as  what  is  seen  when  there  is 
an  increase  in  temperature. 

2.4  Rock  Properties 

Fractured  rock  is  a  very  intricate  system  with  many  fractures  and  complex  matrices  in  which  the 
DNAPL  can  penetrate  deep  into  the  system.  There  are  several  important  properties  that  extend 
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beyond  the  specific  description  of  the  rock  type  such  as  porosity,  compressibility,  thermal 
conductivity,  heat  capacity  and  the  effects  of  temperature  and  pressure  on  these  properties. 

2.4.1  Porosity 

Porosity  is  defined  as  the  volume  of  the  voids  or  open  pore  spaces  divided  by  the  total  volume  of 
the  media  under  consideration.  There  are  two  types  of  porosities:  primary  and  secondary.  The 
primary  porosity  is  created  when  the  rock  is  formed  and  the  secondary  porosity  is  created  due  to 
geological  processes  such  as  weathering,  fracturing,  etc.  (Singhal  and  Gupta,  1999).  The  porosity 
value  is  controlled  by  many  different  factors  such  as  the  “shape  and  arrangement  of  constituent 
grains,  degree  of  sorting,  compaction  and  cementing,  fracturing,  and  solution”  (Singhal  and 
Gupta,  1999).  The  total  porosity  of  the  rock  is  the  sum  of  the  matrix  porosity  and  the  fracture 
porosity  and  is  shown  in  the  equation  below  (Singhal  and  Gupta,  1999): 

0  =  0m  +  0f  (2-8) 

The  porosity  of  the  matrix  (primary  porosity),  0m,  is  the  matrix  void  volume  divided  by  the  total 
volume  where  the  matrix  voids  are  formed  due  to  intergranular  void  spaces  (Singhal  and  Gupta, 
1999).  The  fracture  porosity  (secondary  porosity),  0)  is  equal  to  the  fracture  void  volume  divided 
by  the  total  bulk  volume  (Singhal  and  Gupta,  1999).  The  porosity  can  vary  greatly  for  different 
types  of  rock  and  it  depends  on  the  rock  location  as  well  as  historical  exposure  and  geological 
events. 

2.4.2  Compressibility 

The  total  effective  compressibility  for  fractured  rock  systems  consists  of  three  types  of 
compressibility:  rock  matrix  compressibility,  rock  bulk  compressibility  and  pore  compressibility. 
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All  three  must  be  taken  into  consideration  when  heat  is  being  applied  to  the  system,  as  the 
compressibility  properties  will  be  changed  due  to  changes  in  pressure  (Fatt,  1958).  A  change  in 
the  internal  stress  can  occur  when  there  is  a  decrease  in  the  amount  of  fluid  in  the  pore  spaces  of 
the  rock.  This  stress  change  results  in  changes  in  the  matrix,  bulk,  and  pore  compressibility  in  the 
rock  (Jalalh,  2006).  The  bulk  and  pore  compressibility  are  especially  dependent  on  pressure 
because  if  the  pressure  is  low  then  the  compressibility  will  be  greater  than  if  the  pressure  was 
high  (Geertsma,  1957).  Compressibility  is  also  a  factor  of  the  porosity  as  it  can  vary  greatly 
between  different  rock  types.  When  a  rock  is  subject  to  confining  pressure  and  internal  pore 
pressures,  compression  of  the  rock  occurs. 

If  there  is  a  change  in  pressure  it  results  in  a  fractional  change  of  the  solid  rock  material  (i.e.,  the 
grains),  which  is  known  as  the  rock  matrix  compressibility  (also  known  as  the  grain 
compressibility)  (Amyx  et  al,  1960;  Muggeridge  et  al,  2005;  Geertsma,  1957).  This 
compressibility  is  the  least  important  compared  to  the  rock  bulk  compressibility  and  the  pore 
compressibility.  When  the  porosity  is  greater  than  0. 15,  the  value  of  the  rock  bulk  compressibility 
becomes  much  larger  than  the  matrix  compressibility,  and  thus  it  is  considered  negligible  with 
respect  to  the  bulk  compressibility  (Geertsma,  1957). 

The  rock  bulk  compressibility  occurs  when  there  is  a  “fractional  change  in  volume  of  the  bulk 
volume  of  the  rock  with  a  unit  change  in  pressure”  (Amyx  et  al,  1960,  pg  58).  The  bulk  volume  of 
the  rock  is  the  volume  if  the  pores  are  ignored  (Zimmerman,  1991).  This  compressibility  is 
important  to  consider  if  investigating  under-saturated  reservoirs  because  if  a  pressure  is  applied  to 
the  system  then  the  rock  will  expand  because  it  is  directly  proportional  to  pressure  (Hall,  1953).  A 
large  expansion  of  the  rock  can  alter  the  performance  predictions  of  a  reservoir  (Hall,  1953). 
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The  third  compressibility  to  consider  is  the  pore  compressibility.  When  a  change  in  pressure 
occurs,  there  is  a  change  in  the  pore  volume  of  the  rock  (Hall,  1953).  Laboratory  experiments 
performed  by  Jalalh  (2006)  showed  that  as  the  porosity  of  the  rock  increases,  the  pore  volume 
compressibility  exponentially  decreases.  If  the  porosity  of  the  rock  is  very  low,  then  this  results  in 
higher  pore  compressibility  and  makes  the  variable  important  for  consideration  when  examining 
fractured  rock  systems  (Geertsma,  1957). 


Bulk  compressibility  and  pore  compressibility  of  bedrock  are  described  by  Zimmerman  (1991) 
as: 


£ 


ij 


-1 

V~ 


dV, 


dPi 


]iPk 


(2.9) 


where  s  is  the  compressibility,  i  represents  either  the  bulk  volume  or  pore  volume,  and  j  is  the 
stress  on  the  rock,  which  is  either  confining  pressure  or  the  internal  pore  pressure.  The  bulk 
volume  is  defined  as  the  volume  of  the  pores  and  the  rock  minerals  while  the  pore  volume  is  the 
portion  of  the  bulk  volume  that  is  not  occupied  by  the  rock  minerals.  V  is  the  volume  of  the  rock 
and  dV/dPj  is  the  change  in  volume  due  to  the  changing  pressure  with  the  secondary  pressure,  k, 
being  held  constant. 


There  are  internal  and  external  stresses  present  in  rock  bulk  and  pore  compressibility.  The 
internal  stress  is  usually  hydrostatic  because  it  is  usually  due  to  a  fluid.  The  pore  pressures  are 
varied  and  the  external  pressure  is  held  constant  (Geertsma,  1957).  The  external  stress  can  vary  in 
different  directions  as  it  is  due  to  the  rock  stresses  and  is  when  the  stress  on  the  outer  surface  of 
the  rock  is  changing  and  the  internal  (fluid)  pressure  is  held  constant  (Geertsma,  1957). 
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2.4.3  Thermal  Conductivity 


Thermal  conductivity  is  a  property  that  is  affected  by  the  change  in  several  other  properties  such 
as  porosity,  permeability,  mineralogy  and  specific  surface  area  and  temperature  (Midttoomme 
and  Roaldset,  1998).  It  is  defined  by  Fourier’s  Law,  which  states  that  the  heat  flux,  Q,  is  equal  to 
the  thermal  conductivity,  k  multiplied  by  the  temperature  gradient,  VT  . 

Q  =  -kVT  (2.10) 

Porosity  plays  an  important  role  when  determining  the  thermal  conductivity.  Laboratory  studies 
performed  by  Brigaud  and  Vasseur  (1989),  Thomas  et  al.  (1973)  and  Midttoomme  and  Roaldset 
(1998)  show  that  as  the  porosity  of  a  material  increases,  the  thermal  conductivity  decreases. 

The  decrease  in  thermal  conductivity  is  greater  when  the  porosity  is  smaller.  As  the  porosity 
increases,  the  change  in  thermal  conductivity  decreases,  displaying  a  non-linear  trend  (Brigaud 
and  Vasseur,  1989).  The  thermal  conductivity  of  a  fluid  or  air  is  less  than  that  of  the  grain 
material,  thus  as  the  porosity  increases,  there  are  more  voids  for  the  fluid  or  air  to  reside  in  which 
results  in  a  decrease  of  thermal  conductivity  (Cermack  and  Rybach,  1982;  Midttoomme  and 
Roaldset,  1998). 

Temperature  and  pressure  also  play  an  important  role  on  the  change  in  thermal  conductivity. 
Changes  in  pressure  cause  the  pores  and  microcracks  in  the  rock  to  close  resulting  in  more 
thermal  contact  resistance  (when  grains  come  into  contact  they  create  a  path  of  heat  flow 
resistance),  which  increases  the  thermal  conductivity  (Cermack  and  Rybach,  1982;  Midttoomme 
and  Roaldset,  1998).  The  increase  of  thermal  conductivity  occurs  when  lower  pressures  are 
applied  to  the  rock  because  of  rapid  changes  when  the  cracks  initially  close  (Abdulagatov  et  al, 
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2006).  Once  higher  pressures  are  applied,  the  change  in  thermal  conductivity  is  insignificant 
(Abdulagatov  et  al,  2006). 

The  thermal  conductivity  is  inversely  proportional  to  the  temperature  of  the  rock  (Cermack  and 
Rybach,  1982).  The  increase  in  temperature  causes  a  decrease  in  the  thermal  conductivity.  The 
dependence  on  temperature  is  largely  a  factor  of  the  skeletal  properties  of  the  rock  material  such 
as  the  mineral  components  (Abdulagatov  et  al,  2006). 

2.4.4  Heat  Capacity 

The  specific  heat  capacity,  cp,  is  the  amount  of  heat  that  is  required  to  increase  the  temperature  of 
a  material  by  1  degree  Celsius.  The  heat  capacity  varies  by  temperature  and  pressure.  Equation 
2.7  above  shows  the  equation  for  the  heat  capacity  as  expressed  by  Cermak  and  Rybach  (1982). 


2.5  Thermal  Conductive  Heating 

Thermal  conductive  heating  (TCH),  also  known  as  in  situ  thermal  desorption  (1STD),  is  an 
emerging  technology  in  the  groundwater  remediation  field.  The  remediation  process  uses  a  series 
of  wells  with  heating  elements  inserted  in  each  well.  An  electrical  current  is  provided  at  ground 
surface  to  the  heating  elements.  Figure  2-4  shows  a  fully  installed  and  operable  1STD  setup. 
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Figure  2-4  -  Complete  setup  of  the  thermal  conductive  heating  technology  (TerraTherm, 

Inc.). 

The  heating  element  can  reach  temperatures  of  800  degrees  Celsius.  Resistive  heating  takes  place 
in  the  well  and  the  process  of  conductive  heating  occurs  away  from  the  well  (Baston  and  Kueper, 
2009).  This  process  is  highly  effective  for  chlorinated  solvents  as  the  dominant  removal 
mechanism  is  by  vaporization  (Heron  et  al,  2009).  When  the  chemical  is  heated,  vapors  (steam) 
are  created,  which  are  extracted  by  a  vacuum  system.  Wells  are  typically  placed  6  to  12  ft  apart  to 
ensure  complete  heating  of  the  source  zone  when  dealing  with  semi-volatile  organic  compounds 
and  12  to  20  ft  for  volatile  organic  compounds  (Johnson  et  al,  2009). 

TCH  is  not  as  dependent  on  the  fracture  spacing  for  the  heat  to  be  distributed  as  it  is  in  steam 
injection.  Another  benefit  of  1STD  compared  to  that  of  electrical  resistance  heating  is  that  the 
heating  temperature  can  exceed  100  degrees  Celsius,  which  is  a  major  constraint  of  ERH  (Kunkel 
et  al,  2006).  Also  of  note  is  the  fact  that  ERH  cannot  be  applied  in  many  rock  types  because  of 
the  excessively  high  electrical  resistivity  of  the  rock  matrix. 
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When  TCH  is  applied,  volatile  organic  compounds  are  destroyed  by  the  following  mechanisms  as 
described  by  Baker  et  al.  (2006):  (i)  evaporation,  (ii)  steam  distillation,  (iii)  boiling,  (iv) 
hydrolysis,  (v)  oxidation,  and  (vi)  pyrolysis. 

With  any  remediation  technology,  groundwater  sources  can  pose  a  problem  as  an  inflow  of 
groundwater  can  reduce  the  heating  in  the  system.  To  reduce  the  heat  losses,  the  incoming 
groundwater  can  be  boiled  by  either  steam  injection  or  placement  of  an  extra  row  of  heater  wells 
outside  of  the  treatment  zone  to  preheat  the  water  (Baston  and  Kueper,  2009).  To  reduce  the 
amount  of  cooling  that  occurs,  especially  when  high  temperatures  are  desired,  an  insulation 
blanket  will  be  placed  over  the  surface  of  the  treatment  area  in  the  form  of  a  concrete  slab  (Baston 
and  Kueper,  2009).  This  not  only  reduces  heat  loss,  but  also  reduces  the  amount  of  infiltration  of 
precipitation  into  the  system. 

Thermal  conductive  heating  has  been  successfully  implemented  at  contaminated  sites.  A  study 
done  by  Triplet  Kingston  et.  al  (2010)  found  that  from  1988  to  2007  there  have  been  26  pilot/field 
scale  TCH  applications  with  17  of  these  being  completed  since  2000.  The  study  concluded  that 
there  has  been  only  one  field  scale  application  performed  in  weathered  bedrock  with  the 
remaining  applications  conducted  in  porous  media  settings. 

Examples  of  TCH  applications  include: 

•  Baker  et  al,  (2006)  discussed  a  former  manufactured  gas  plant  site  in  North  Adams, 
Massachusetts  in  which  residual  coal  tar  was  remediated  with  TCH  inside  of  a  brick 
walled  gasholder  of  5.5  m  depth  which  was  backfilled  with  a  mixture  of  silt,  sand  gravel 
and  cobbles  and  debris.  The  contaminates  observed  include  naphthalene,  benzo(a)pyrene, 
benzene  and  other  Petroleum  hydrocarbons  (PHCs).  The  gasholder  was  heated  for  3 
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months  at  low  temperatures  and  then  temperatures  of  325°C  were  sustained  for  9  months. 
There  was  a  96.43  to  99.98%  reduction  in  concentration  measured  in  the  treatment  zone. 

•  Heron  et  al.  (2008)  used  TCH  to  remediate  a  site  in  southeastern  U.S.  which  had  a 
trichloroethylene  source  zone  that  extended  90  ft  below  the  ground  surface.  The 
treatment  zone  geology  consisted  of  25  ft  of  fill  underlain  by  30  ft  of  severely  weathered 
granitic  gneiss,  20  ft  of  weathered  bedrock  and  fractured  bedrock  at  the  bottom  of  the 
treatment  zone.  The  treatment  zone  was  heated  at  an  average  temperature  of  100°C  for 
110  days.  Starting  soil  and  rock  concentrations  were  measured  to  be  as  high  as 
81,000,000  pg/kg  and  was  reduced  to  17  pg/kg  post  treatment.  The  groundwater 
concentration  in  the  treatment  zone  was  reduced  by  74.5  to  99.7%. 

•  Baker  et  al.  (2007)  discussed  the  remediation  of  a  contaminated  site  in  Southern 
California  Edison,  U.S.  The  treatment  zone  consisted  primarily  of  silty  soils  which  were 
contaminated  with  polycyclic  aromatic  hydrocarbons  (PAHs),  pentachlorophenol  (PCP), 
and  polychlorinated  dibenzodioxins  and  furans  (PCDD/Fs).  The  treatment  zone  was  a 
maximum  depth  of  32  m  and  the  soil  concentrations  for  PAHs,  PCP  and  PCDD/Fs  were 
measured  to  be  an  average  of  30,600,  2,940  and  18  pg/kg  pre -treatment  respectively.  The 
subsurface  was  heated  for  30  days  at  300°C  which  resulted  in  reduced  soil  concentrations 
to  59,  1,250  and  0.1 1  pg/kg  for  PAHs,  PCP  and  PCDD/Fs  respectively  which  was  below 
the  remedial  goal  for  the  site. 
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Chapter  3 


Numerical  Simulation  of  Trichloroethylene  Remediation  in  Fractured 
Bedrock  by  Thermal  Conductive  Heating 


Abstract 

A  thermal  conductive  heating  (TCH)  pilot  test  was  conducted  in  2009  at  the  Naval  Air  Warfare 
Center  (NAWC)  in  West  Trenton,  New  Jersey  which  was  historically  impacted  by 
trichloroethylene  (TCE).  Rock  samples  were  collected  pre-  and  post-remediation  inside  of  the  test 
pilot  area  to  assess  TCE  removal  from  the  rock  matrix.  Using  site-specific  data  collected  during 
the  pilot  test,  a  two-dimensional  finite  difference  model  was  created  using  TMVOC  to  investigate 
TCE  removal  from  the  rock  matrix  compared  to  the  removal  during  the  NAWC  test  pilot.  In  the 
baseline  model  simulation,  TCE  removal  in  the  rock  matrix  simulation  was  67%  compared  to 
63.5%  removal  in  the  pilot  test.  The  numerical  model  was  used  as  a  screening  tool  to  examine 
what  parameters  affect  the  mass  removal  from  the  rock  matrix  during  1 00  days  of  heating.  A 
limited  sensitivity  analysis  examined  the  matrix  permeability,  matrix  porosity  and  energy  rate  to 
the  subsurface.  The  sensitivity  analyses  indicated  that  TCE  mass  removal  increases  with 
increased  matrix  porosity  and  increased  rate  of  energy  application. 

Keywords:  Thermal  Conductive  Heating,  fractured  bedrock,  trichloroethylene, 
remediation,  numerical  model 
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3.1  Introduction 


The  remediation  of  dense,  non-aqueous  phase  liquids  (DNAPLs)  such  as  creosote,  coal  tar  and 
chlorinated  solvents  in  bedrock  is  a  challenging  task  given  the  discrete  nature  of  flow  through 
fractures  and  the  potential  for  back  diffusion  of  contaminants  in  many  rock  types.  Technologies 
that  rely  upon  fluid  injection,  such  as  steam  injection,  in-situ  chemical  oxidation,  enhanced  in-situ 
bioremediation,  and  co-solvent  flushing  are  typically  limited  by  sweep  efficiency  and  a  general 
inability  to  address  contaminants  sequestered  in  the  low  permeability  porous  rock  matrix  in  a 
timely  manner.  Thermal  conductive  heating  (TCH),  which  relies  largely  upon  heat  conduction 
rather  than  the  delivery  of  injected  fluids,  has  potential  applicability  in  fractured  rock  because  it 
can  deliver  heat  to  not  only  open  fractures,  but  to  the  rock  matrix  as  well. 

TCH  employs  electrical  heaters  suspended  inside  a  cased  borehole  to  deliver  energy  to  the 
surrounding  formation  and  can  achieve  temperatures  of  up  to  approximately  500  °C  (Heron  et  al., 
2009).  The  heat  migrates  away  from  the  heater  borings  by  a  combination  of  thermal  conduction 
(driven  by  a  temperature  gradient)  and  convection  (migration  of  steam  produced  by  boiling 
ground  water).  Upon  heating  to  the  boiling  point  of  water,  steam  is  formed  in-situ  along  with  the 
transfer  of  contaminant  mass  to  the  vapor  phase  (Heron  et  al.,  2009).  Vapor  and  liquid  recovery 
is  typically  achieved  through  vacuum  extraction  and  conventional  groundwater  pumping, 
respectively.  In  addition  to  physical  removal  through  vapor  and  liquid  recovery,  biological  and 
chemical  degradation  mechanisms  may  occur  during  and  after  thermal  remediation.  These 
mechanisms  may  include  thermal  destruction  by  oxidation  and  pyrolysis  near  heating  elements  at 
temperatures  of  approximately  400°C,  microbial  mineralization  of  NAPL  components,  and 
hydrolysis  at  elevated  temperatures  (Baker  and  Kuhlman,  2002).  The  required  target  temperature 
in  the  treatment  zone  as  well  as  the  treatment  time  varies  depending  on  a  number  of  factors 

including  the  properties  of  the  target  compounds  and  heat  losses.  Heat  losses  generally  occur  at 
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the  surface,  below  and  around  the  perimeter  of  the  treatment  zone  and  are  influenced  by 
groundwater  flux  (Baston  and  Kueper,  2009).  To  reduce  heat  loss,  the  heater  borings  are 
generally  extended  1  to  2  m  outside  of  the  treatment  zone  (Johnson  et  al.,  2009).  The  general 
physical  layout  of  a  TCH  application,  including  depiction  of  a  power  supply,  produced  fluids 
treatment  components,  temperature  monitoring  points,  and  pressure  monitoring  points  is 
illustrated  in  Figure  3-1. 


Granular  Activated 


Figure  3-1  -  Typical  configuration  for  thermal  conductive  heating  (Johnson  et  al.,  2009, 

Appendix  A:  Heron  and  Baker). 


TCH  has  been  implemented  on  many  occasions  in  porous  media  (Triplett  Kingston  et  al.,  2010), 
but  there  is  relatively  little  experience  in  fractured  rock.  To  evaluate  the  performance  of  TCH  in 
fractured  rock,  a  field  pilot  test  was  carried  out  by  TerraTherm,  Inc.  in  2009  in  collaboration  with 
the  Naval  Facilities  Engineering  Services  Center,  Queen’s  University  and  the  United  States 
Geological  Survey  at  the  former  Naval  Air  Warfare  Center  (NAWC)  in  West  Trenton,  New 
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Jersey.  The  NAWC  facility  operated  from  the  mid-1950s  to  the  late  1990s  where  it  served  as  a 
U.S.  Navy  jet  engine  testing  facility.  Activities  on  the  site  resulted  in  the  release  of 
trichloroethylene  (TCE)  to  bedrock. 

Using  data  collected  at  the  site  during  the  pilot  test  and  associated  site  characterization  activities, 
a  two-dimensional  (2D)  radial  finite  difference  model  was  developed  using  TMVOC  (Pruess  and 
Battistelli,  2002).  TMVOC  is  part  of  the  TOUGH2  family  of  codes  and  is  capable  of  simulating 
multiphase  flow,  heat  transfer,  and  transport  of  volatile  organic  compounds  (VOCs)  in  three- 
dimensional  heterogeneous  porous  media  or  fractured  rock  (Pruess  and  Battistelli,  2002).  The  2D 
model  was  developed  as  a  screening  tool  to  evaluate  the  reduction  of  TCE  in  the  model  to  that  of 
the  pilot  test  when  similar  heating  conditions  were  applied.  A  limited  sensitivity  analysis  was 
conducted  to  focus  on  how  the  performance  of  TCE  extraction  was  influenced  by  changes  in 
matrix  permeability,  matrix  porosity  and  energy  input  into  the  subsurface.  The  results  from  this 
study  will  generate  recommendations  for  future  TCH  remediation  applications  in  fractured 
bedrock  aquifers. 

3.2  TCH  Pilot  Test 

A  TCH  pilot  test  was  performed  by  TerraTherm  Inc.  at  the  NAWC  facility  between  April  10  and 
July  24,  2009.  A  network  of  15  heater/extraction  wells  (vapor  and  liquid  extraction)  was  installed 
over  an  area  of  45.6  nr.  The  wells  were  installed  to  a  depth  of  16.8  m  below  ground  surface 
(bgs)  and  screened  their  entire  length  to  allow  for  maximum  extraction  of  fluids.  Further  details 
pertaining  to  the  pilot  test  are  described  in  Rodriguez  et  al.  (2012).  Figure  3-2  illustrates  the 
configuration  of  heater/extraction  wells  and  temperature  monitoring  points  associated  with  the 
field  pilot  test. 
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Figure  3-2  -  Plan  view  configuration  of  heater/extraction  wells  and  temperature  monitoring 

points. 


Three  wells  were  installed  on  the  upgradient  side  of  the  heating  network  to  aid  in  reducing  heat 
losses  due  to  incoming  groundwater  flow.  A  total  of  eight  temperature  monitoring  borings  (T1  to 
T8)  were  completed  with  10  temperature  measurement  points  in  each  spaced  at  1.52  m  intervals 
beginning  at  1.52  m  below  ground  surface.  A  concrete  cap  was  placed  on  top  of  the  treatment 
zone  to  reduce  heat  losses  at  ground  surface.  The  pilot  test  included  combined  heating  and  fluid 
extraction  for  98  days,  followed  by  fluid  extraction  alone  for  4  days  post-heating.  Steady-state 
temperatures  achieved  in  the  treatment  zone  varied  from  approximately  60  °C  to  approximately 
110  °C  with  an  average  of  98  °C  (Figure  3-3).  Temperatures  in  the  bottom  4.6  m  of  the  test 
interval  did  not  reach  100  °C,  likely  because  of  the  cooling  influence  of  inflowing  groundwater  at 
those  depths. 
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Figure  3-3  -  Average  temperature  for  eight  temperature  monitor  borings  (T1-T8)  that  had 
1.52  m  interval  temperature  points  from  1.52  m  to  15.24  m  below  ground  surface  during  the 
NAWC  pilot  test  (data  provided  by  TerraTherm,  Inc.). 

To  evaluate  the  removal  of  TCE  from  the  rock  matrix,  three  boreholes  (BR-1,  BR-2,  and  BR-3) 
were  drilled  in  the  centre  of  the  treatment  area  in  which  samples  were  collected  pre-treatment 
(Figure  3-2).  Three  additional  boreholes  were  drilled  and  sampled  post-treatment  at  locations 
adjacent  to  the  pre-treatment  boreholes.  Rock  samples  were  collected  every  0.3  to  1.5  m  with 
BR-1  containing  the  greatest  density  of  sampling  pre-  and  post-treatment.  Rock  samples  were 
placed  in  a  rock  crushing  device,  preserved  in  methanol  and  shaken  for  24  hours  to  allow  the  TCE 
to  partition  into  the  methanol  (Rodriguez  et  al.,  2012).  The  methanol  solution  was  then  extracted 
from  each  rock  sample  and  analyzed  for  TCE  concentration  by  an  external  laboratory.  Figure  3  A 
presents  the  rock  matrix  sampling  results  as  a  function  of  depth  from  the  BR-1  pair  of  coring 
locations  pre-  and  post-treatment.  This  location  exhibited  higher  initial  concentrations  of  TCE  in 
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the  rock  matrix  than  BR-2  and  BR-3.  Considering  results  from  all  three  paired  locations, 
approximately  63.5%  of  the  TCE  was  removed  from  the  rock  matrix  as  a  result  of  TCH 
application. 


0  50  100  150  200  250  300 

TCE  (mg/kg) 

Figure  3-4  -  TCE  concentration  in  the  rock  matrix  from  location  BR-1  pre-  and  post¬ 
treatment. 


3.3  Model  Development 
3.3.1  Site  Characterization  Input 

The  NAWC  Facility  is  underlain  by  mudstone  bedrock  that  is  part  of  the  Lockatong  Formation 
(Tiedeman  et  al.,  2010;  Lacombe  and  Burton,  2010).  Figure  3-5  below  displays  a  cross  section  of 
the  geology  at  the  NAWC  Facility  (Tiedeman  et  al.,  2010)  with  the  pilot  test  area  highlighted. 
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Figure  3-5  -  Cross  Section  of  the  geology  at  the  NAWC  site  with  the  pilot  test  location 

shown  (Tiedeman  et  al.,  2010). 

The  overburden  material  consists  of  a  thin  layer  of  fill  and  weathered  saprolite  which  varies  from 
1  to  8  m  in  thickness  across  the  NAWC  site.  In  the  area  of  the  pilot  test,  the  depth  of  this  unit  was 
to  approximately  1.8  m  below  ground  surface  (bgs).  The  saprolite  is  composed  primarily  of  silt 
and  clay  and  has  a  bulk  average  permeability  of  6.5x  1  O' 13  m2  (Tiedman  et  al.,  2010). 

This  layer  is  underlain  by  a  highly  fractured,  weathered  grey  mudstone  in  which  “the  individual 

mudstone  beds  are  difficult  to  distinguish”  (Tiedeman  et  al.,  2010).  The  weathered  mudstone  is 

underlain  by  dipping  formations  that  alternate  between  beds  of  highly  fractured  mudstones  which 

are  highly  conductive  due  to  the  inclined  fractures;  and  more  competent,  less  fractured  beds 

(Tiedeman  et  al.,  2010;  Lewis-Brown  et  al.,  2006).  Laboratory  testing  completed  on  the  mudstone 

rock  in  the  area  of  the  pilot  test  (from  approximately  9. 1  m  bgs)  was  measured  to  have  a  bulk 

density  of  2520  kg/m3,  a  matrix  porosity  of  3.3%,  a  matrix  permeability  of  1x10  17  m2  and  a 

matrix  fraction  organic  carbon  of  0.0079  (average  properties  from  Rodriguez  et  al.,  2012). 

Hydraulic  testing  was  completed  by  TerraTherm  Inc.  prior  to  the  pilot  test  heating.  This  was  done 
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by  employing  a  straddle  packer  assembly  and  was  conducted  in  the  three  holes  in  which  the  rock 
samples  were  taken.  The  bulk  permeability  of  the  pilot  test  area  was  measured  to  be  8.4x10  14  m2. 
The  bulk  permeability  is  three  orders  of  magnitude  greater  than  the  matrix  permeability  as  the 
bulk  permeability  takes  into  consideration  possible  fracture  networks  in  the  weathered  grey 
mudstone  whereas  the  laboratory  testing  only  considers  a  small  rock  core  disk. 

The  water  table  across  the  NAWC  Facility  varies  from  1  to  5  m  below  ground  surface  (Tiedeman 
et  al.,  2010)  but  in  the  area  of  the  TCH  pilot  test,  the  water  table  was  located  1.52  m  bgs.  The 
hydraulic  gradient  at  the  site  in  the  pilot  test  area  was  approximately  0.008  (Tiedman  et  al.,  2010). 

3.3.2  Model  Domain  and  Boundary  Conditions 

TMVOC  was  used  to  create  a  two-dimensional  (2D)  radial  finite  difference  model  to  simulate 
TCE  removal  from  the  rock  matrix  using  data  collected  from  the  pilot  test.  Figures  3-6  and  3-7 
display  the  conceptual  2D  model  upon  which  the  numerical  model  is  based.  The  numerical 
simulations  were  designed  to  be  radially  symmetric  to  simulate  the  heater/extraction  well  that  is 
located  in  the  center  of  the  well  network  (Figure  3-6). 

The  inside  radius  of  the  model  domain  is  a  symmetry  boundary.  The  column  of  cells  directly 
adjacent  to  this  boundary  (the  first  column  of  cells)  has  heating/extraction  properties  assigned. 
The  outer  radius  of  the  domain  was  placed  1.5  m  away  and  is  subjected  to  no-flow  conditions  in 
order  to  model  the  symmetry  boundary  of  a  repeated  pattern  which  represents  the  first  ring  of 
wells  in  the  field  test  pilot. 

The  observation  location  is  the  outermost  column  of  cells,  which  is  located  1 .25  m  away  (at  the 
centroid  of  the  cell  column)  from  the  center  well.  This  location  in  the  model  represents  the 
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location  half  way  between  the  inner  well  and  first  ring  of  wells  in  the  field  pilot  test  which  is 
approximately  where  the  pre-  and  post-remediation  rock  samples  were  collected  (BR-1  to  BR-3). 


•  Heating  and  Extraction  Well 
O  Observation  Location 
-  Model  Domain 
r  Radius  =  1.5  m 


Figure  3-6  -  Plan  view  of  the  conceptual  model. 
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Figure  3-7  -  Cross-sectional  view  of  the  conceptual  2D  radial  symmetric  model. 
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The  top  boundary  was  fixed  with  a  constant  temperature  of  10  °C  and  a  constant  pressure  of 
101,325  Pa  representing  atmospheric  conditions  while  the  bottom  boundary  was  set  to  a  no-flow 
boundary.  The  bottom  boundary  was  extended  10.5  m  below  the  treatment  zone  in  order  to 
reduce  the  influence  of  the  bottom  boundary  condition  on  the  model  results.  Before  setting  the 
bottom  boundary  to  a  no-flow  boundary,  a  constant  pressure  of  356,231  Pa  was  applied  and  the 
model  was  run  to  steady  state  in  order  to  create  a  water  table  approximately  1.5  m  below  the  top 
of  the  model.  One  limitation  of  the  radially  symmetric  domain  is  that  a  hydraulic  gradient  cannot 
be  incorporated  into  the  simulations.  This  presents  ideal  conditions  as  there  was  a  hydraulic 
gradient  observed  during  the  pilot  test. 

Three  horizontal  fractures  were  incorporated  into  the  model  domain  as  illustrated  in  Figure  3-7. 
As  discussed  in  Section  3.3.1,  the  bedrock  underlying  the  NAWC  site  is  composed  of  dipping 
layers  of  highly  weathered  and  variably  fractured  rock.  Dipping  fractures  could  not  be 
incorporated  into  the  model  domain  due  to  the  finite  element  grid  selection  as  well  as  being 
radially  symmetrical.  The  three  fracture  locations  were  selected  based  on  the  sharp  spikes  in  TCE 
concentration  found  in  the  rock  matrix  pre -remediation  samples  as  shown  in  Figure  3-4.  The 
fracture  zone  cells  had  a  1  cm  thickness  and  were  assigned  a  permeability  calculated  according  to 
the  following  depth-weighted  arithmetic  mean  (Baston  et  al.,  2010): 


kfz 


kf(e)  +  km(zfz-e) 


" fz 


(3.1) 


where  kfz  is  the  permeability  assigned  to  the  fracture  zone  cells,  kf  is  the  fracture  permeability 
calculated  as  e2/l 2,  e  is  the  fracture  hydraulic  aperture  assumed  to  be  equal  to  the  mechanical 
aperture,  km  is  the  matrix  permeability,  and  zf-  is  the  thickness  of  the  fracture  zone  cell  (1  cm). 
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The  fracture  hydraulic  aperture  was  calculated  on  the  basis  of  hydraulic  testing  completed  at  the 
site  by  TerraTherm,  Inc.  in  borings  BR-1  through  BR-3.  A  total  of  32  hydraulic  aperture  values 
were  calculated  ranging  from  54  microns  to  234  microns  with  an  arithmetic  mean  value  of  138 
microns.  The  average  aperture  was  used  in  Equation  3.1  to  calculate  the  permeability  of  the 
fracture  zone  cells  and  assigned  to  each  of  the  three  fractures  in  the  model  domain.  A  matrix 
permeability  of  l.OxlO'17  m2  was  assigned  to  the  mudstone,  which  represents  the  arithmetic  mean 
of  two  measurements  performed  on  rock  disks  cut  from  pre -remediation  rock  core  obtained  from 
within  the  treatment  zone.  The  resulting  fracture  zone  permeability  is  calculated  to  be  l.OxlO'11 
m2. 

The  TCE  profile  used  as  the  baseline  concentration  for  the  model  simulations  was  based  on  the 
average  of  the  pre -remediation  rock  matrix  concentrations  measured  during  the  pilot  test  in  BR- 1 , 
BR-2  and  BR-3. 

3.3.3  Baseline  Model  Establishment 

The  primary  objective  of  the  numerical  simulations  was  to  compare  rock  matrix  TCE  mass 
removal  achieved  during  the  pilot  test  to  that  predicted  by  the  model  when  subjected  to  similar 
heating  conditions  (i.e.  temperature  history).  The  average  temperature  recorded  at  measurement 
points  T1  through  T8  during  the  pilot  test  were  matched  by  adjusting  the  daily  energy  input  into 
the  model.  The  matrix  permeability  was  initially  set  to  lxl O'17  m2  and  the  fracture  zone 
permeability  set  to  1x10  11  m2  to  match  the  average  temperature  profile. 

Once  the  temperature  profile  was  calibrated,  the  percentage  of  mass  removal  from  the  rock  matrix 
was  estimated  to  be  54%  which  is  less  than  the  average  TCE  removal  measured  during  the  pilot 
test  at  63.5%.  The  pilot  test  at  the  NAWC  site  was  conducted  in  a  highly  fractured,  weathered 
grey  mudstone  which  would  allow  for  liquid  and  vapour  extraction  through  the  fracture  network. 
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The  model  domain  incorporated  three  fractures  without  vertical  connections,  as  it  assumes 
competent  bedrock  matrix  in  between  the  fractures.  To  better  simulate  the  subsurface  conditions 
within  the  pilot  test  zone  (that  likely  contained  several  interconnected  micro-fractures  creating 
additional  pathways  for  vapour  extraction),  a  sensitivity  analysis  was  completed  using  a  range  of 
matrix  permeability  values  from  1x10  13  m2  to  lxlO'17  m2.  When  changing  the  matrix  permeability 
from  lxlO'13  m2  to  lxlO'17  m2  the  average  temperature  in  the  model  at  the  observation  location 
varied  by  approximately  4  degrees.  The  greater  the  matrix  permeability  the  lower  the  average 
temperature  profile.  This  is  shown  in  Figure  3-8.  The  permeability  of  the  fracture  zone  cells  was 
not  recalculated  using  the  equations  above  when  the  model  matrix  permeability  was  varied. 


Figure  3-8  -  Temperature  versus  time  curve  displaying  the  influence  of  matrix  permeability. 
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The  resulting  percentage  of  TCE  removed  from  the  rock  matrix  during  the  various  permeability 
simulations  is  displayed  Table  3-1.  This  table  shows  that  the  percentage  of  TCE  removed  is 
greatest  when  the  matrix  permeability  is  set  to  1x10  15  m2,  which  had  a  67%  removal  as  compared 
to  the  pilot  test  of  63.5%.  The  model  bulk  permeability  was  calculated  using  the  depth  weighted 
arithmetic  mean  of  the  model  matrix  permeability  and  the  model  fracture  permeability.  When  the 
model  matrix  permeability  was  set  to  lxlO'15  m2,  the  model  bulk  permeability  was  calculated  to 
be  3.5x10  14  m2  which  is  comparible  to  the  bulk  permeability  measured  during  hydraulic  testing 
completed  by  TerraTherm  Inc.  in  2009  (a  value  of  8.4xl0'14  m2  was  estimated). 


Table  3-1  -  Matrix  bulk  permeability  of  the  numerical  simulations  with  different  matrix 

permeabilities  applied. 


Model  Matrix 
Permeability 
(m2) 

Model  Fracture 
Permeability 
(m2) 

Model  Bulk 
Permeability 
(m2) 

Percent  TCE 
Removed 
(%) 

lxlO'17 

l.OxlO'11 

2.0xl0"14 

54 

lxlO'16 

l.OxlO"11 

2. lxlO"14 

54 

lxlO'15 

l.OxlO'11 

3.5xl0'14 

67 

lxlO'14 

l.OxlO"11 

1.7xl0"13 

39 

lxlO'13 

l.OxlO'11 

1.5xl0"12 
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Intuitively,  the  greatest  amount  of  extraction  was  anticipated  to  occur  with  the  highest  matrix 
permeability  of  lxlO'13  m2,  but  this  was  not  the  case.  The  primary  route  of  TCE  extraction  occurs 
in  the  three  fractures.  Closer  examination  of  the  TCE  concentrations  at  different  time  steps  was 
completed  with  the  matrix  permeability  set  to  lxlO"15  m2  and  1  x  1 O'1’  m2. 


Figure  3-9  displays  the  model-predicted  TCE  profile  in  the  rock  matrix  at  the  observation  location 

over  time  with  a  matrix  permeability  of  lxlO'15  m2.  This  figure  shows  that  as  the  temperature  in 

the  rock  matrix  increases,  boiling  begins  to  occur  in  the  fractures.  At  30  days  the  matrix 

temperature  was  an  average  of  80  °C  (prior  to  boiling)  and  the  TCE  concentration  in  the  fracture 
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zone  increased  by  2  to  8  times  greater  than  the  initial  concentration.  This  increase  is  a  result  of 
the  TCE  partitioning  to  the  vapour  phase  in  the  small  fracture  zone  areas.  This  phenomenon  has 
been  observed  by  Chen  et  al.  (2010)  during  rock  core  laboratory  experiments  when  boiling 
occurred  in  the  rock  core. 
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Figure  3-9  -  Model-predicted  TCE  profile  in  the  rock  matrix  at  the  observation  location 
over  time  when  the  matrix  permeability  is  set  to  lxlO"15  m2.  Grey  curves  represent  the  initial 

TCE  concentration  specified  in  the  model. 
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The  TCE  concentration  reduces  primarily  in  the  fracture  zone  cells  of  the  upper  two  fractures 
after  heating  for  40  days  due  to  extraction  of  the  TCE  vapours  and  liquids.  Once  heating  of  the 
matrix  has  occurred  for  50  days,  extraction  is  taking  place  in  all  three  fractures.  The  concentration 
above  the  upper  fracture  begins  to  increase  in  concentration.  By  80  days  the  TCE  concentration 
between  the  top  and  bottom  fracture  is  significantly  decreased,  but  the  concentration  above  the 
top  fracture  has  increased  in  concentration.  Once  1 00  days  of  heating  has  occurred  (the  end  of  the 
numerical  simulation),  the  TCE  profile  shows  that  the  TCE  concentration  above  the  upper 
fracture  has  increased.  The  TCE  below  the  bottom  fracture  shifted  its  location  slightly  downward 
with  the  final  concentration  remaining  nearly  the  same  as  the  initial  concentration  in  the  bottom 
fracture.  This  shift  in  concentration  was  observed  when  the  matrix  permeability  was  set  to  lxlO'13 
nr  as  discussed  below. 

Figure  3-10  shows  the  model-predicted  TCE  profile  in  the  rock  matrix,  at  the  observation  location 
over  time,  when  the  matrix  permeability  is  set  to  1x10  13  nr.  Boiling  begins  to  occur  in  the  model 
after  30  days  of  heating.  Because  the  matrix  permeability  is  only  2  orders  of  magnitude  smaller 
than  the  fracture  zone  permeability  (1x10  11  m2),  the  TCE  partitions  into  the  vapour  phase  not 
only  in  the  fracture  zone  cells,  but  also  in  the  matrix.  The  TCE  mobilized  differently  than 
previously  observed.  The  concentration  increased  in  the  area  above  the  uppermost  fracture  and  in 
between  the  uppermost  and  middle  fractures. 

At  40  days,  boiling  continues  and  the  TCE  vapour  in  between  the  uppermost  fracture  and  bottom 
fracture  are  extracted,  therefore  creating  a  significant  concentration  decrease  in  the  matrix.  TCE 
vapour  migrates  to  above  and  below  the  fracture  zone  creating  an  increase  in  concentration.  As 
the  simulation  continues  to  1 00  days,  the  concentration  above  the  top  fracture  decreased  slightly, 
but  the  TCE  continues  to  migrate  deeper  into  the  rock  matrix  below  the  bottom  fracture. 
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Depth  (m) 


Figures  3-9  and  3-10  demonstrate  that  as  the  rock  matrix  permeability  is  increased,  becoming 
closer  in  magnitude  to  the  value  of  the  fracture  zone  permeability,  the  TCE  is  able  to  mobilize 
beyond  the  fracture  zone  cells.  This  occurs  due  to  the  pressure  gradients  that  are  formed  when  the 
permeability  is  altered.  Figure  3-1 1  displays  the  pressure  profiles  at  the  observation  location  for 
the  two  different  scenarios  discussed  above. 
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Matrix  Permeability  =  IE-13  m2 
Matrix  Permeability  =  IE-15  m2 
Fracture 
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Figure  3-11  -  Model-predicted  pressure  profiles  at  the  observation  location  when  the  matrix 
permeability  is  set  to  1x10'°  m2  and  lxlO'15  m2. 
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As  previously  discussed,  the  primary  route  of  TCE  extraction  is  through  the  fracture  zone  cells. 
When  the  temperature  in  the  rock  matrix  is  increased,  the  pressure  is  also  increased.  For  all  time 
steps  shown  in  Figure  3-1 1  when  the  matrix  permeability  is  set  to  lxlO'13  nr  (therefore  closer  in 
magnitude  to  the  fracture  zone  cells),  the  pressure  distribution  is  relatively  constant  between  the 
fracture  zone  cells  and  the  matrix.  There  is  a  pressure  decrease  towards  the  fill  and  saprolite  layer 
as  well  as  towards  the  bottom  of  the  heating  zone. 

When  the  permeability  is  decreased  (1x10  13  nr),  pressure  gradients  are  formed  towards  the 
fracture  zone  areas  except  for  above  the  top  fracture  where  the  pressure  gradient  is  driven  towards 
the  higher  permeability  fill  area. 

These  pressure  gradients  dictate  the  redistribution  of  TCE  seen  in  Figures  3-9  and  3-10.  When  the 
TCE  is  vapourized  it  is  driven  by  the  pressure  gradients  and  will  migrate  from  a  high  pressure 
area  to  a  low  pressure  area.  When  the  matrix  permeability  is  lxlO"13  nr,  there  is  no  significant 
pressure  gradient  towards  the  fracture  zone  cells  but  instead  above  and  below  the  fracture  area 
and  thus  there  was  redistribution  of  TCE  above  and  below  the  fractures.  When  the  permeability  is 
set  to  lxl0'1?  m2,  there  are  4  orders  of  magnitude  different  between  the  fracture  zone  cells  and  the 
matrix  therefore  the  pressure  gradients  are  created  and  the  TCE  vapours  travel  from  the  higher 
pressures  in  the  matrix  to  the  lower  pressures  in  the  fracture  zone  where  it  is  ultimately  removed. 
Because  of  the  higher  permeability  fill  layer  above  the  fracture  area,  there  is  a  gradient  created 
towards  the  fill  and  hence  some  TCE  migration  to  the  zone  above  the  uppermost  fracture. 

From  this  sensitivity  analysis,  it  can  be  concluded  that  the  optimal  matrix  permeability  was  1x10" 
15  m2  as  it  has  the  greatest  percentage  of  TCE  removed  from  the  rock  matrix  and  did  not  cause 
significant  redistribution  of  the  TCE  above  and  below  the  heating  zone.  For  these  reasons,  a 
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matrix  permeability  of  lxl O' 15  m2  was  used  for  the  baseline  model  simulation  to  complete  the 
sensitivity  analysis  discussed  in  Section  3.3.4. 


The  model  simulations  represent  an  idealized  scenario  by  having  three  primary  fractures. 
Increasing  the  matrix  permeability  helps  to  account  for  the  possible  fracture  networks  in  the 
subsurface  at  the  NAWC  site.  The  simulations  clearly  demonstrate  that  the  primary  pathway  for 
TCE  removal  is  through  the  fracture  zone  cells  and  when  the  matrix  permeability  becomes  close 
to  the  fracture  zone  permeability  the  driving  force  (the  pressure)  towards  the  fractures  is 
significantly  decreased. 

A  summary  of  the  baseline  model  input  parameters  used  for  the  sensitivity  analysis  discussed  in 
Section  3.3  is  provided  in  Table  3-2.  Thermal  conductivity,  specific  heat  capacity,  pore 
compressibility  and  the  Klinkenberg  parameter  were  assumed  to  be  equal  for  the  saprolite  and 
bedrock.  Scaled  power  functions  were  used  by  TMVOC  to  calculate  the  gas  and  liquid  water 
phase  relative  permeabilities  following  Stone  (1970).  The  capillary  pressure  functions  were 
adapted  from  Falta  et  al.  (1992),  which  were  calculated  using  modified  versions  of  the  van 
Genuchten  (1980)  constitutive  model  given  by  Parker  et  al.  (1987).  The  total  mass  of  TCE 
removed  from  the  rock  matrix  during  the  pilot  test  was  estimated  to  be  63.5%  of  the  initial 
amount  (Rodriguez  et  al.,  2012).  The  calibrated  baseline  model  predicted  67%  TCE  mass 
removal. 
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Table  3-2  -  Baseline  model  input  parameters. 


Parameter 

Units 

Saprolite 

Mudstone 

Bedrock 

Fracture  Zone 

Cells 

Density 

kg/m3 

2520 

2520 

2520 

Permeability 

m2 

6.5xl0'13 

l.OxlO'15 

l.OxlO'11 

Porosity1 

% 

30 

3.3 

3.3 

C  2 

'-'wet 

W/m°C 

2.84 

2.84 

2.84 

c  2 

'-'dry 

W/m°C 

2.19 

2.19 

2.19 

Specific  heat  capacity3 

J/kg°C 

1000 

1000 

1000 

Pore  compressibility4 

Pa1 

3.0xl0'9 

3. 0x1  O'9 

3. 0x1  O'9 

Klinkenberg  parameter5 

Pa1 

2.4x1  O'6 

2.4x1  O'6 

2.4x1  O'6 

Fraction  organic  carbon 

- 

0.0079 

0.0079 

0.00 

1.  The  saprolite  layers  were  assumed  to  have  a  porosity  of  30%  which  is  a  common  porosity  for  silts  and 
clays  as  specified  by  Hough  (1957). 


2.  Cwet  is  the  thermal  conductivity  under  saturated  conditions  and  is  the  thermal  conductivity  under 
de-saturated  conditions.  This  parameter  was  calculated  using  the  equations  in  Woodside  and  Messmer 
(1961)  with  a  thermal  conductivity  of  the  matrix  of  3.0  W/m°C. 

3.  Scharli  and  Rybach  (2001);  Cermak  and  Rybach  (1982). 

4.  Zimmerman  (1991). 

5.  Webb  and  Pruess  (2003). 


3.3.4  Sensitivity  Analysis  Outline 

The  baseline  calibrated  model  was  used  to  perform  a  sensitivity  analysis  to  examine  what 
parameters  are  most  influential  in  dictating  the  amount  of  TCE  mass  removal  from  the  rock 
matrix.  Table  3-3  summarizes  the  parameters  that  were  varied  in  addition  to  the  matrix 
permeability  in  the  sensitivity  analysis  and  the  range  of  values  that  were  considered. 
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Table  3-3  -  Sensitivity  analysis  parameters. 


Parameter 

Field  Parameters 

Unit 

Baseline  Model 

Values  Considered 

Permeability 

m2 

lxlO'15 

10'17,  1016,  10'14,  10'13 

Matrix  Porosity 

% 

3.3 

1,5,  10 

Thermal  Conductive  Heating 
Parameter 

Energy  Rate 

J/s 

Field  Temperature 
Match 

-10%,  +10%,  +20%,  +30% 

3.4  Sensitivity  Results  and  Discussion 
3.4.1  Matrix  Porosity 

Primary  porosity  is  created  when  the  rock  is  formed  and  is  referred  to  here  as  the  matrix  porosity; 
secondary  porosity  is  created  due  to  geological  processes  such  as  weathering  and  fracturing 
(Singhal  and  Gupta,  1999).  The  2D  simulations  only  considered  primary  porosity  with  the 
baseline  model  having  a  matrix  porosity  of  3.3%,  which  was  determined  through  laboratory 
testing  on  core  samples  obtained  from  within  the  target  heating  zone. 

Figure  3-11  displays  the  model-predicted  temperature -time  curves  when  the  matrix  porosity  is 
altered  between  1  and  10%.  Figure  3-12  displays  the  percentage  of  TCE  removed  as  a  function  of 
matrix  porosity.  Figures  3-11  and  3-12  illustrate  that  the  matrix  porosity  does  not  significantly 
affect  the  temperature  of  the  rock  matrix  when  the  same  amount  of  energy  is  applied  in  the 
numerical  model.  Flowever,  the  matrix  porosity  does  influence  the  percentage  of  TCE  removed 
from  the  rock  matrix.  As  the  porosity  increases,  the  TCE  removed  at  100  days  also  increases. 
There  was  10%  more  TCE  removed  from  the  rock  matrix  when  the  matrix  porosity  was  increased 
from  3.3%  to  10%. 
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Figure  3-13  -  Percentage  of  TCE  removed  from  the  rock  matrix  at  100  days  when  the 
matrix  porosity  is  varied  from  1  to  10%. 
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3.4.2  Energy  Rate 


A  variable  energy  rate  (J/s)  was  applied  to  the  column  of  matrix  blocks  where  the 
heater/extraction  well  was  located  on  the  inner  radius  of  the  model  domain.  A  higher  energy  was 
applied  until  the  maximum  average  temperature  was  reached  in  which  the  energy  was  then 
reduced  to  sustain  that  temperature.  The  baseline  model  had  an  energy  input  that  created  a 
temperature-time  profile  similar  to  what  was  recorded  during  the  NAWC  pilot  test.  This  energy 
was  then  increased  or  decreased  to  examine  what  effects  the  energy  input  into  the  rock  matrix  has 
on  TCE  removal. 

The  baseline  model  energy  input  was  varied  from  a  10%  decrease  in  energy  to  a  30%  increase  in 
energy,  as  shown  in  Figure  3-13.  When  the  temperature  in  the  rock  matrix  is  increased  by 
increasing  the  energy  input,  the  percentage  of  TCE  removed  is  also  increased  (Figure  3-14).  A 
decrease  in  energy  rate  by  10%  reduces  the  temperature  to  near  the  boiling  temperature  of  TCE 
and  this  resulted  in  20%  of  TCE  being  removed  from  the  rock  matrix. 

When  a  30%  increase  in  energy  input  is  applied,  the  maximum  temperature  achieved  in  the  rock 
matrix  was  112  °C  and  86%  of  the  TCE  removed.  While,  when  a  20%  increase  in  energy  was 
applied  it  resulted  in  84.5%  TCE  removal.  It  can  be  concluded  from  this  study  that  a  greater 
energy  into  the  rock  matrix  results  in  higher  maximum  temperatures  reached  in  the  model  but 
there  is  little  final  benefit  to  increasing  the  energy  input  to  more  than  20%  as  the  final 
contaminant  removal  after  100  days  between  20%  and  30%  was  insignificant. 
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Figure  3-14  -  Temperature  versus  time  curve  displaying  the  influence  of  a  change  in  energy 

input  to  the  rock  matrix. 


Figure  3-15  -  Percentage  of  TCE  removed  from  the  rock  matrix  at  100  days  when  the 
energy  input  into  the  2D  simulations  is  varied. 
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3.5  Conclusions 


A  thermal  conductive  heating  field  pilot  test  was  conducted  in  2009  at  the  Naval  Air  Warfare 
Center  in  West  Trenton,  New  Jersey  in  which  historical  TCE  contamination  is  found  in  fractured 
mudstone  bedrock.  The  pilot  test  consisted  of  a  network  of  15  heater/extraction  wells  that  were 
operated  for  102  days  to  remove  TCE  from  the  rock  matrix.  Using  the  data  collected  at  the  site,  a 
2D  screening  model  was  created  using  TMVOC  in  which  the  primary  objective  was  to  examine 
the  TCE  mass  removal  in  the  rock  matrix  from  the  numerical  simulations  compared  to  the  pilot 
test. 

The  pilot  test  was  found  to  have  a  TCE  removal  in  the  rock  matrix  of  63.5%  as  compared  to  the 
baseline  model  simulation  which  had  a  removal  of  67%.  The  numerical  model  presented  ideal 
conditions,  such  as  only  3  fractures  separated  by  competent  bedrock  and  no  hydraulic  gradient.  A 
sensitivity  analysis  on  the  rock  matrix  permeability  was  conducted  to  calibrate  the  model  with  the 
field  results,  ft  was  concluded  that  for  the  numerical  simulations  a  matrix  permeability  of  lxlO"1''’ 
m2  was  appropriate  due  to  the  idealized  model  conditions.  The  fracture  zone  cells  were  identified 
to  be  the  primary  area  where  TCE  removal  occurred  due  to  the  higher  permeability  in  these  cells 
compared  to  the  lower  permeability  in  the  rock  matrix.  The  sensitivity  analysis  concluded  that  the 
greater  the  maximum  temperature  sustained  in  the  model,  the  greater  the  percentage  of  TCE 
removed  from  the  rock  matrix. 

Future  investigations  using  numerical  simulations  could  create  a  model  with  a  more  realistic 
fracture  network  in  the  subsurface  rather  than  ideal  conditions.  Creating  a  3 -dimensional  model 
rather  than  a  radially  symmetric  model  would  also  allow  for  a  hydraulic  gradient  to  be  applied  to 
the  model.  Using  a  2D  model  does  not  account  for  the  heating  influence  of  nearby  wells, 

therefore  creating  a  3D  model  would  also  allow  for  these  interactions  to  be  observed.  Heating 
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well  interactions  could  also  increase  the  subsurface  temperatures  therefore  resulting  in  greater 
TCE  removal  from  the  rock  matrix. 
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Chapter  4 


Conclusions  and  Recommendations 

In  2009  a  pilot  test  was  conducted  at  the  Naval  Air  Warfare  Center  in  West  Trenton,  New  Jersey 
which  employed  thermal  conductive  heating.  The  site  is  historically  contaminated  with 
trichloroethylene  in  the  mudstone  bedrock.  The  pilot  test  was  completed  on  a  small  portion  of  the 
site  with  an  area  of  45.6  nr.  A  network  of  15  heater/extraction  wells  that  conducted  vapour  and 
liquid  extraction  was  installed  to  a  depth  of  16.8  m  below  the  ground  surface.  Each  well  was 
screened  to  allow  for  maximum  extraction  of  fluids  and  the  pilot  test  was  conducted  over  102 
days. 

Using  the  data  collected  during  the  pilot  test,  a  two-dimensional  numerical  model  was  created 
using  TMVOC  which  is  part  of  the  TOUGH2  family  of  codes  that  is  capable  of  simulating 
multiphase  flow,  heat  transfer,  and  transport  of  volatile  organic  compounds.  The  primary 
objective  was  to  create  a  screening  model  in  order  to  examine  the  TCE  mass  removal  from  the 
rock  matrix  in  the  numerical  summations  compared  to  the  pilot  test. 

Three  rock  core  samples  were  taken  pre-  and  post-remediation  during  the  pilot  test  at  the  NAWC 
site.  The  core  samples  were  measured  for  TCE  concentrations  before  and  after  heating  to  measure 
the  percentage  of  mass  removed  from  the  rock  matrix.  The  pilot  test  resulted  in  an  average  of 
63.5%  TCE  removal  after  102  days.  The  two-dimensional  radial  finite  difference  model 
simulations  were  designed  to  simulate  the  heater/extraction  well  that  was  located  in  the  center  of 
the  well  network  during  the  pilot  test.  The  simulations  contained  three  primary  fracture  zone 
locations  which  were  selected  based  upon  the  pre -remediation  rock  core  TCE  concentration 
findings  with  competent  bedrock  in  between.  Temperature  data  collected  during  the  pilot  test  was 
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used  in  the  numerical  simulations  to  match  the  heating  patterns  observed.  The  baseline  model  had 
a  67%  TCE  removal  from  the  rock  matrix  after  100  days  of  heating  which  was  similar  to  that 
measured  during  the  pilot  test. 

The  numerical  model  presents  ideal  conditions  as  it  only  contained  three  fractures  separated  by 
competent  bedrock.  As  it  was  a  radially  symmetric  model  no  hydraulic  gradient  can  be 
incoiporated.  The  pilot  test  was  conducted  in  highly  fractured  weathered  bedrock  which  most 
likely  contained  numerous  fracture  networks.  To  help  account  for  this,  the  matrix  permeability  in 
the  numerical  simulations  was  increased.  The  simulations  demonstrated  that  the  primary  pathway 
for  extraction  was  through  the  three  fracture  zones  in  the  model. 

A  limited  sensitivity  analysis  was  completed  on  the  baseline  simulation  which  also  examine  the 
effect  that  matrix  porosity  and  energy  input  into  the  subsurface  has  on  the  effectiveness  of  the 
remedial  system.  The  sensitivity  analysis  concluded  that  the  greater  the  maximum  temperature 
sustained  in  the  model,  the  greater  the  percentage  of  TCE  removed  from  the  rock  matrix.  It  also 
demonstrated  that  there  was  not  a  significant  benefit  (i.e.  TCE  removal)  when  a  30%  increase  in 
energy  was  applied  to  the  model  versus  20%  for  100  days.  The  results  from  this  study  may  help 
generate  recommendations  for  future  TCH  remediation  applications  in  fractured  bedrock  aquifers 
as  the  simulations  demonstrated  that  the  fracture  zone  plays  an  important  role  in  the  amount  of 
contaminant  removed  from  the  rock  matrix. 

Future  investigations  using  TMVOC  could  be  completed  to  further  examine  the  success  of  mass 
removal  from  fractured  rock.  Future  investigations  could  include: 

•  Incorporating  a  fracture  network  in  the  subsurface  rather  than  three  primary  fractures 
which  account  for  vertical  and  horizontal  connections. 


66 


•  Creating  a  three-dimensional  model  rather  than  a  radially  symmetric  model  would  allow 
for  a  hydraulic  gradient  to  be  applied  to  the  model  which  may  also  affect  the  amount  of 
contaminant  removal  from  the  subsurface. 

•  Incorporation  the  15  heating/extraction  wells  into  a  3D  model  domain  in  a  similar  layout 
as  the  pilot  test.  This  would  allow  for  the  influence  and  interactions  of  the  nearby  wells  in 
the  well  network  which  were  not  accounted  for  in  the  two-dimensional  model  to  be 
examined. 

•  Conduct  a  sensitivity  analysis  on  the  horizontal  and  vertical  grid  discritization. 

•  Conduct  a  sensitivity  analysis  on  the  permeability  of  the  fracture  zone  cells  to  examine  if 
there  is  more  TCE  extraction  with  increased  fracture  zone  permeability.  This 
permeability  could  be  varied  using  the  range  of  fracture  apertures  calculated  from  the 
hydraulic  testing  during  the  pilot  test. 
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Appendix  A 


TMVOC  Governing  Equations 

This  section  provides  a  brief  summary  of  the  governing  equations  and  assumptions  that  are  used 
by  TMVOC,  A  Numerical  Simulator  for  Three-Phase  Non-isothermal  Flows  in  Multicomponent 
Hydrocarbon  Mixtures  in  Saturated-Unsaturated  Heterogeneous  Media.  Further  information  can 
be  found  in  the  TMVOC  user  manual  by  Pruess  and  Battistelli  (2002). 

TMVOC  is  a  numerical  model  developed  at  the  Lawrence  Berkeley  National  Laboratory,  which 
is  capable  of  simulating  three-dimensional  heterogeneous  porous  media  or  fractured  rock.  It  is  an 
extension  of  the  TOUGH2  simulator  and  considers  multi-component  volatile  organic  compound 
(VOC)  flow  in  all  three  phases  (gas-water-NAPL)  above  and  below  the  groundwater  table. 
TMVOC  is  able  to  simulate  non-condensable  gases  as  well  as  advection,  dispersion,  soiption, 
phase  partitioning  and  biodegradation.  TMVOC  allows  for  irregular  grids  to  be  created  by  using 
the  Integral  Finite  Difference  Method  (1FDM)  for  space  discritization. 

The  structure  and  execution  of  TMVOC  uses  the  same  architecture  as  TOUGH2  (Figure  A-l). 
Primary  and  secondary  input  variables  are  specified  by  the  user  for  which  the  numerical  model 
uses  a  series  of  linear  equations  and  equations  of  state  (EOS)  modules  to  produce  the  output. 
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■EOS-Modute' 


Figure  A-l  -  TMVOC  architecture  (Pruess  and  Battistelli,  2002). 


TMVOC  can  simulate  heat  transfer  applications  such  as  thermal  remediation  and  steam  injection 
remediation  technologies.  Heat  transfer  simulated  by  TMVOC  occurs  by  conduction  and  multi¬ 
phase  convection.  TMVOC  assumes  that  the  three  phases  are  in  “local  chemical  and  thermal 
equilibrium”. 


Air  and  other  non-condensable  gases  can  be  incorporated  into  the  numerical  simulator.  If  there 
are  VOC  vapors  present  in  the  system  then  Raoult’s  Law  governs  which  is  dependent  on  the 
aqueous  phase  concentration  and  vapor  pressure  for  single  component  mixtures,  which  can  be 
seen  in  Equation  A.  1 . 


Xg  Pg  —  xn^vap 


(A-l) 
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where  x  is  the  mole  fraction  of  the  chemical  component  in  the  gas  (g)  and  the  NAPL  (//),  and  P  is 
the  pressure  of  the  gas  and  the  vapor  pressure. 

The  aqueous  phase  is  dominated  by  water.  When  VOCs  are  present,  the  composition  of  the  water 
is  dependent  on  the  concentration  in  NAPL  or  gas  phase  and  the  solubility  of  the  VOC.  When 
non-condensable  gases  are  present  they  are  represented  by  TMVOC  using  Henry’s  Law. 


The  NAPL  phase  (oil  phase)  considers  only  VOCs  and  assumes  that  the  solubility  of  water  and 
non-condensable  gases  (NCGs)  into  the  NAPL  phase  is  minimal.  The  solubility  of  the  NAPL  in 
water  follows  the  equation  given  by  Hoot  et  al.  (1957): 


x%  =  0.281226  -  4.23610U  +  1.75244u2 


(A.2) 


where, 


v  =  log10(1.87c  +  32) 


(A- 3) 


and  Tc  is  the  temperature  in  degrees  Celsius. 

TMVOC  uses  a  general  mass  and  energy  balance  equation  as  follows: 


d  f 

Tt  J  "rdli 


I  Ft  ■  ndT„  +  j 


bydvn 


(A.4) 


where  M  is  the  mass  accumulation  term  for  the  component  y  (if  f3  1  it  denotes  water,  NCGs  or 
VOCs;  or  y=  NK  +  1  to  represent  the  heat  component),  and  is  integrated  over  a  domain  V„.  It  is 
equal  to  the  integral  of  the  mass  or  heat  flux,  F,  “bounded  by  the  closed  surface”,  Tn,  plus  the 
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integration  of  the  sinks  and  sources,  b.  “n  is  a  normal  vector  on  surface  element  dfn,  pointing 


inward  into  Vn”. 


The  mass  accumulation  term  is  represented  for  the  water  and  NCG  components  in  the  model  by 
Equation  A. 5. 


MY  =  r\  ^  SpPpXYp 
P 


(A.  5) 


where  y  can  be  the  mass  or  heat  component  for  the  fluid  phases  (3,  which  can  be  in  a  liquid,  gas  or 
NAPL  state,  q  is  the  porosity,  S  is  the  saturation,  p  is  the  density  and  x  is  the  mole  fraction.  The 
mass  and  heat  accumulation  formulas  include  rock  properties  such  as  sorption  onto  the  rock 
grains  (represented  in  second  half  of  Equation  A.6)  and  the  specific  heat  of  rock,  CR  for  the  heat 
accumulation  equation  (Equation  A.  7). 


MK  =  cj)^  SpppXp  +  (1  -  < p)pRpwxYwKd 
P 


(A.6) 


with  p  representing  the  density  (of  the  phase,  (3,  rock,  R,  or  water,  w)  and  Kd  the  distribution 
coefficient.  The  heat  accumulation  formula  is  as  follows: 


MNK+ 1  =  (1  -  0)pRCRT  +  0  £  Spppup 

P 


(A.7) 


where  T  is  the  temperature  and  Up  is  the  specific  internal  energy.  The  accumulation  terms  are 
discretized  in  space  using  the  1FDM. 
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Appendix  B 


Model  Testing 

When  using  an  existing  numerical  model,  it  is  important  to  determine  if  it  has  been  correctly 
compiled.  To  determine  that  TMVOC  was  compiled  correctly,  Example  Problem  5  from  the 
document  entitled  “TMVOC,  A  Numerical  Simulator  for  Three-Phase  Non-isothermal  Flows  of 
Multicomponent  Hydrocarbon  Mixtures  in  Saturated-Unsaturated  Heterogeneous  Media”  written 
by  Pruess  and  Battistelli  (2002)  was  used  and  is  further  referred  to  as  the  TMVOC  Users  Manual. 
The  TMVOC  Users  Manual  provides  a  copy  of  the  input  files  which  can  be  used  directly  to  run 
the  simulations.  This  example  was  selected  because  it  not  only  contained  input  files  which  could 
be  used  to  run  the  simulations,  but  several  graphed  figures  that  could  be  re-generated  for 
comparison  to  verify  the  model  was  compiled  correctly. 

Example  Problem  5  simulates  a  laboratory  experiment  conducted  by  Hunt  et  al.  (1988)  which  was 
“designed  to  evaluate  the  mobilization  and  transport  of  a  [non  aqueous  phase  liquid  (NAPL)] 
during  water  and  steam  flooding”  (Falta  et  al.,  1992).  The  simulation  modeled  a  91  cm  horizontal 
laboratory  column  with  a  cross  sectional  area  of  20.428  cnr  and  examines  the  displacement  of  a 
benzene-toluene  mixture  during  waterflooding.  The  model  consisted  of  50  evenly  space  elements 
filled  with  sand  which  had  a  density  of  2650  kg/m3,  a  porosity  of  38.5%  and  a  permeability  of 
1.6x10  11  m2.  The  first  part  of  the  simulation  was  to  inject  the  benzene-toluene  mixture  into  the 
15th  element  (approximately  27  cm  from  the  inlet)  at  a  constant  rate  into  the  column  over  a  120 
second  period.  Figure  B-l  shows  the  conceptual  model  of  the  laboratory  experiment  and  model 
simulations  (Falta  et  al,  1992). 
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Injection 

(inlet) 


NAPL  Injection  Point 


Production 

(outlet) 


Figure  B-l  -  Conceptual  diagram  of  the  model  simulations  (Falta  et  al.,  1992). 


The  second  step  of  the  simulation  was  to  induce  a  waterflood  into  the  horizontal  column. 
Waterflooding  is  a  process  by  which  water  at  22  °C  is  injected  into  the  left  side  of  the  horizontal 
column  and  allowed  to  flow  through  the  column  and  drain  out  the  opposite  side.  The  simulation 
injected  10.2  pore  volumes  (PV)  of  water  into  the  column  for  the  waterflood. 


Figure  B-2  below  shows  the  non  aqueous  phase  liquid  (NAPL)  saturations  along  the  91  cm 
horizontal  column  distance.  The  figure  displays  the  NAPL  saturation  before  waterflood  which  is 
part  1  of  the  simulation  when  the  benzene-toluene  mixture  was  injected  for  120  s.  The  waterflood 
occurred  for  the  second  part  of  the  simulation  and  the  NAPL  saturation  of  the  benzene-toluene 
mixture  was  shown  for  2,422  seconds  which  corresponds  to  1.2  pore  volumes  (PV),  18,162 
seconds  (9  PV)  and  108,  971  s  (10.2  PV). 
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Figure  B-2  -  Re-created  Figure  10.5.5  from  the  TMVOC  Users  Manual,  Example  Problem 

No.  5  to  demonstrate  model  calibration. 


To  further  examine  if  calibration  of  TMVOC  was  correct,  Figure  10.5.6  from  Example  Problem  5 
in  the  TMVOC  Users  Manual  was  recreated  and  is  shown  in  Figure  B-3.  Figure  B-3  displays  the 
benzene  and  toluene  aqueous  mass  fractions  over  time  at  the  horizontal  column  outlet  (at  91  cm). 
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Figure  B-3  -  Re-created  Figure  10.5.6  from  the  TMVOC  Users  Manual,  Example  Problem 

No.  5  to  demonstrate  model  calibration. 


Figure  B-3  shows  that  the  amount  of  benzene  removed  at  the  end  of  the  waterflood  is  greater  than 
the  mass  fraction  of  toluene  removed.  This  occurs  because  the  solubility  of  benzene  is  greater 
than  toluene  which  allows  the  benzene  to  have  more  mobility. 

Recreation  of  the  Figures  10.5.5  and  10.5.6  from  Example  Problem  5  in  the  TMVOC  users 
manual  was  successful  in  testing  that  the  model  was  compiled  correctly  (as  shown  by  the 
matching  Figures  B-2  and  B-3). 
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Appendix  C 


TMVOC  Model  Input  Property  Definitions 

C.l  Rock  Properties 

Density  [kg/m3]  -  The  mass  of  a  substance  per  unit  volume  of  the  substance. 

Permeability  [m2]  -  Property  of  the  soil  or  rock  that  dictates  the  ability  of  a  fluid  to  move  by 
advection.  The  permeability  depends  greatly  on  the  pore  size  and  pore  connectivity. 

Fraction  of  Organic  Carbon  -  The  fraction  of  organic  carbon  that  is  found  naturally  occurring 
in  the  subsurface.  It  is  measured  as  mass  of  organic  carbon  per  mass  of  soil  or  rock. 

Formation  Heat  Conductivity  [W/m-C]  -  Also  known  as  the  thermal  conductivity  of  a  material, 
it  is  quite  variable  for  different  rock  types.  It  is  highly  dependent  on  saturation  and  temperature 
changes.  The  heat  conductivity  can  be  for  liquid  saturated  or  desaturated  conditions.  The  reader  is 
referred  to  Chapter  2.0  Literature  Review  for  further  details. 

Rock  Grain  Specific  Heat  [J/kg-C]  -  Also  known  as  the  specific  heat  capacity,  is  the  amount  of 
energy  that  is  required  to  increase  the  temperature  of  a  material  by  1  degree  Celsius. 

Porosity  -  The  ratio  of  volume  of  the  void  space  of  a  porous  material  to  the  total  volume  of  the 
material.  Changes  in  porosity  occur  when  there  is  a  change  in  temperature  and  pressure.  It  is 
defined  by  Reid  et  al,  (1987)  and  TMVOC  uses  the  following  equation  (Falta  et  al.,  1992): 
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0  =  0r[1+sp(P-Pr)  +  et(T-Tr)] 


(C.1) 


where  0R  is  the  reference  porosity  of  the  material,  8P  is  the  pore  compressibility  and  ST  is  the 

expansivity  and  are  both  assumed  to  be  constant.  P  and  T  represent  the  pressure  and  temperature 
being  applied  to  the  system  with  PR  and  TR  being  the  reference  pressure  and  temperature 
respectively. 

Pore  Compressibility  [Pa'1]  -  When  there  is  change  in  pressure,  compression  of  the  pores  in  the 
rock  matrix  occurs.  The  reader  is  referred  to  Chapter  2.0  Literature  Review  for  further  details.  It 
is  implemented  into  the  model  using  the  following  (Pruess  &  Battistelli,  2002): 


£p 


(C.2) 


where  8P  is  the  pore  compressibility,  0  is  the  porosity,  and  30/ dP  is  the  change  in  porosity  with 
respect  to  pressure  at  a  constant  time,  T. 


Pore  Expansivity  [°C  ']  -  The  change  in  porosity  when  expansion  occurs  in  response  to  a  change 
in  temperature  with  a  constant  pressure.  TMVOC  applies  the  following  equation  (Pruess  & 
Battistelli,  2002): 


sT  — 


(C.3) 


where  sT  is  the  pore  expansivity,  and  30 /dT  is  the  change  in  porosity  with  respect  to  temperature 
at  a  constant  pressure,  P. 
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Tortuosity  Factor  -  The  actual  travel  distance  through  a  porous  medium  divided  by  the  straight- 
line  travel  distance.  When  the  tortuosity  factor  for  binary  diffusion  is  set  to  0,  TMVOC  calculates 


the  gas  phase  tortuosity,  Tg,  using  the  Millington  and  Quirk  Model  (1961)  which  uses  the 
porosity,  0,  and  the  gas  phase  saturation,  Sg  as  follows  (Falta  et  al.,  1992): 


zg  =  0'/V/3 


(C.4) 


Klinkenberg  Parameter  -  The  slip  flow  of  a  gas  in  a  pore,  which  reduces  the  drag  forces  and  in 
turn  enhances  the  gas  flow  through  the  pores  (Webb,  2001).  It  is  important  when  the  permeability 
of  the  media  is  less  than  10  18  m2  and  at  low  differential  pore  pressures  (Webb  and  Preuss,  2003). 
TMVOC  calculates  the  absolute  permeability  of  the  gas  phase  according  to  (Pruess  &  Battistelli, 


2002): 


k  =  kOQ 


(C.5) 


where  kx  is  the  permeability  at  “infinite”  pressure,  b  is  the  Klinkenberg  parameter  and  P  is  the 
pressure.  The  Klinkenberg  parameter  can  be  estimated  using  the  Jones  and  Owens  (1980) 
equation  for  low  permeability  sands  in  the  range  of  10'14  to  10  19  m2  (Webb  &  Pruess,  2003)  as 
follows: 

b  =  k0~° 33  (C.6) 


where  ka  is  the  intrinsic  permeability  in  m2. 
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C.2  Chemical  Properties 


Diffusion  -  Movement  (non-advective)  of  a  solute  from  a  region  of  high  concentration  to  a 
region  of  low  concentration.  The  diffusive  flux  (one-dimensional)  is  represented  by  Fick’s  Law 
(Fetter,  2001): 


F  = 


dC 

dx 


(C.7) 


where  F  is  the  mass  flux  of  the  solute,  D  is  the  diffusion  coefficient  of  the  solute,  C  is  the  solute 
concentration  and  dC/dx  is  the  concentration  gradient. 
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Appendix  D 


Calculations 

D.l  Fracture  Aperture  and  Permeability 

Hydraulic  tests  employing  a  straddle  packer  assembly  were  conducted  in  the  three  well  holes  in 
which  the  rock  samples  were  taken,  BR-1,  BR-2  and  BR-3.  The  hydraulic  testing  was  completed 
by  TerraTherm  Inc.  prior  to  the  pilot  test  heating.  From  the  hydraulic  testing,  TerraTherm  Inc. 
calculated  the  hydraulic  conductivity.  The  fracture  aperture  to  be  incorporated  into  the  TMVOC 
simulation  was  calculated  using  the  results  of  these  hydraulic  tests.  The  measured  transmissivity 
(TR)  of  the  tested  interval  length  (/)  is  related  to  the  bulk  rock  hydraulic  conductivity  ( K)  through 
(Fetter,  2001): 


TR  =  Kl 


(D.l) 


If  it  is  assumed  that  the  rock  matrix  permeability  is  negligible  compared  to  that  of  the  fracture, 
then  the  measured  transmissivity  can  be  related  to  the  fracture  hydraulic  aperture  ( e )  through  (de 
Marsily,  1986): 


TR  =  —  (e)3  (D'2) 

12 iiy  J 

where  p  is  the  water  density,  g  is  gravitational  acceleration,  p  is  the  viscosity  of  water  and  e  is  the 
fracture  aperture.  Rearranging  this  equation,  the  fracture  aperture  can  be  calculated  as  (with 
results  displayed  in  Tables  D-l  to  D-3): 
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e  = 


3)1277?^ 

J  Pd 


(D.3) 


It  is  assumed  here  that  the  groundwater  temperature  is  approximately  1 0  degrees  Celsius  and  thus 
the  density  of  water  was  1  g/cm3  and  the  viscosity  was  0.01  g/s-cm  (Fetter  2001). 


The  rock  permeability  can  also  be  calculated  using  the  packer  testing  results.  The  permeability  is 
calculated  as: 


K P  (D.4) 

P9 
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Table  D-l-  Hydraulic  test  and  aperture  calculation  results  for  BR-1  (hydraulic  conductivity  calculated  by  TerraTherm  Inc.; 
transmissivity  and  fracture  aperture  calculated  by  Ashley  McKenzie). 


Date 

Packer 

Depth 

Well 

Depth 

Packer 

Interval 

b 

Hydraulic 

Conductivity 

K 

Transmissivity 

T 

Fracture 

Aperture 

2b 

Average 

Aperture 

2b 

Permeability 

(ft) 

(ft) 

(ft) 

(cm/s) 

(cm'Vs) 

(cm) 

(pm) 

K) 

3.15E-05 

9.60E-03 

1.05E-02 

3.21E-14 

11/21/2008 

11 

21 

10 

2.29E-05 

6.97E-03 

9.48E-03 

99 

2.33E-14 

2.41  E-05 

7.35E-03 

9.65E-03 

2.46E-14 

2.73E-04 

8.33E-02 

2.17E-02 

2.79E-13 

11/21/2008 

21 

31 

10 

1.96E-04 

5.99E-02 

1.94E-02 

198 

2.00E-13 

1.67E-04 

5.09E-02 

1.84E-02 

1.70E-13 

1.12E-04 

3.42E-02 

1.61  E-02 

1.14E-13 

11/21/2008 

31 

41 

10 

9.90E-05 

3.02E-02 

1.55E-02 

160 

1.01E-13 

1.16E-04 

3.55E-02 

1.63E-02 

1.19E-13 

9.49E-05 

2.89E-02 

1.52E-02 

9.68E-14 

11/21/2008 

41 

51 

10 

7.92E-05 

2.41  E-02 

1.43E-02 

141 

8.07E-14 

5.66E-05 

1.72E-02 

1.28E-02 

5.77E-14 
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Table  D-2-  Packer  test  and  aperture  calculation  results  for  BR-2  (hydraulic  conductivity  calculated  by  TerraTherm  Inc.;  transmissivity 
and  fracture  aperture  calculated  by  Ashley  McKenzie). 


Date 

Packer 

Depth 

(ft) 

Well 

Depth 

(ft) 

Packer 

Interval 

b 

(ft) 

Hydraulic 

Conductivity 

K 

(cm/s) 

Transmissivity 

T 

(cmz/s) 

Fracture 

Aperture 

2b 

(cm) 

Average 

Aperture 

2b 

(pm) 

Permeability 

(m2) 

11/20/2008 

16 

26 

10 

1.97E-04 

5.99E-02 

1.94E-02 

194 

2.01E-13 

11/20/2008 

26 

36 

10 

0.00E+00 

- 

- 

145 

- 

0.00E+00 

- 

- 

- 

5.65E-05 

1.72E-02 

1.28E-02 

5.76E-14 

1.14E-04 

3.47E-02 

1.62E-02 

1.16E-13 

11/20/2008 

36 

46 

10 

0.00E+00 

- 

- 

117 

- 

1.55E-05 

4.74E-03 

8.34E-03 

1.58E-14 

9.24E-05 

2.82E-02 

1.51  E-02 

9.42E-14 

8.82E-06 

2.69E-03 

6.90E-03 

9.00E-15 

11/20/2008 

46 

56 

10 

2.01  E-05 

6.12E-03 

9.08E-03 

85 

2.05E-14 

2.21  E-05 

6.73E-03 

9.38E-03 

2.25E-14 
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Table  D-3  -  Packer  test  and  aperture  calculation  results  for  BR-3  (hydraulic  conductivity  calculated  by  TerraTherm  Inc.;  transmissivity 
and  fracture  aperture  calculated  by  Ashley  McKenzie). 


Date 

Packer 

Depth 

Well 

Depth 

Packer 

Interval 

b 

Hydraulic 

Conductivity 

K 

Transmissivity 

T 

Fracture 

Aperture 

2b 

Average 

Aperture 

2b 

Permeability 

(ft) 

(ft) 

(ft) 

(cm/s) 

(cm'Vs) 

(cm) 

(pm) 

K) 

4.72E-05 

1 .44E-02 

1.21  E-02 

4.82E-14 

11/19/2008 

11 

21 

10 

3.05E-05 

9.29E-03 

1.04E-02 

116 

3.11E-14 

5.03E-05 

1.53E-02 

1.23E-02 

5.13E-14 

0.00E+00 

- 

- 

- 

0.00E+00 

- 

- 

- 

11/19/2008 

21 

31 

10 

5.16E-06 

1.57E-03 

5.77E-03 

83 

5.26E-15 

3.75E-06 

1.14E-03 

5.19E-03 

3.82E-15 

7.05E-05 

2.15E-02 

1.38E-02 

7.19E-14 

2.99E-04 

9.12E-02 

2.23E-02 

3.05E-13 

11/19/2008 

31 

41 

10 

1.65E-04 

5.03E-02 

1.83E-02 

191 

1.68E-13 

1.22E-04 

3.71  E-02 

1.66E-02 

1.24E-13 

6.19E-05 

1.89E-02 

1.32E-02 

6.31E-14 

11/19/2008 

41 

51 

10 

2.58E-05 

7.87E-03 

9.88E-03 

102 

2.63E-14 

1.63E-05 

4.97E-03 

8.47E-03 

1.66E-14 

1.99E-05 

6.06E-03 

9.05E-03 

2.03E-14 
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D.2  Klinkenberg  Coefficient 


The  Klinkenberg  coefficient  was  calculated  using  Equation  D.4  as  given  by  Jones  and  Owens  (1980)  for 
low  permeability  sands  having  a  permeability  value  between  10' 19  and  10  14  m2  (Webb  and  Pruess,  2003). 
The  average  bulk  rock  permeability  (k0)  at  the  Naval  Air  Warfare  Centre  test  site  in  Trenton,  New  Jersey 
is  9.79E-18  nr  (Rodriguez  et  al.,  2012)  and  bair  is  the  Klinkenberg  parameter.  The  resulting  Klinkenberg 
coefficient  used  for  the  2D  numerical  model  simulations  was  2.4  x  10'6  Pa"1. 

Klinkenberg  Coefficient  =  - —  =  ko~033 

bair 


D.3  Rock  Pore  Compressibility 

The  pilot  test  occurred  primarily  in  a  highly  weathered  grey  mudstone  which  is  classified  as  a 
sedimentary  rock.  The  compressibility  of  rock  can  be  highly  variable  and  was  not  measured  during  the 
NAWC  pilot  test.  The  pore  compressibility  of  the  mudstone  was  calculated  using  the  information  from 
Zimmerman  (1991)  on  pages  29-30  (compressibility  values  for  several  sandstone  reservoirs).  The  average 
pore  compressibility  was  calculated  from  twenty-four  reservoir  sandstones.  The  rock  pore 
compressibility,  sp,  is  calculated  using  Equation  D.5. 

_ Cr  ( D  .5 ) 


where  sb  is  the  bulk  compressibility,  Cr  is  the  addition  factor  given  by  Fatt  (1958a)  as  1.8x10  7  psi'1  and  0 
is  the  porosity  of  the  sandstone.  The  average  rock  pore  compressibility  was  calculated  to  be  3xl0"9  Pa"1. 


D.4  Fracture  Zone  Cell  Permeability 


The  permeability  of  the  fracture  zone  cell  can  be  calculated  using  the  permeability  found  from  the  packer 
tests  and  analytical  tests  using: 


T  k^e  (D.6) 

(Azm  +  e) 


where  kb  is  the  bulk  permeability  from  the  packer  testing  which  is  equal  to  8.4x10  14  m2,  km  is  the  matrix 
permeability  from  analytical  testing  set  to  9.79xl0‘18  m2  and  A zm  is  the  packer  spacing  equal  to  3.048  m. 
The  permeability  in  the  field  of  the  fracture  (kj)  is  calculated  using  Equation  D.7  below  and  substituted  in 
Equation  D.6. 


(D.7) 


The  fracture  zone  permeability  (kjz)  used  for  the  numerical  model  simulations  is  calculated  using  Equation 
D.8  as  developed  by  Baston  et  al.  (2010). 


kfz 


kf  '  (?  T  km(Azfz  (?) 


A  z, 


■fz 


(D.8) 


where  A Zyz  is  the  cell  thickness  in  the  model  for  the  fracture  zone  which  was  set  to  1  cm.  The  resulting 
permeability  values  used  in  the  numerical  2D  model  is  lxl 0"17  m2  for  used  as  the  bulk  rock  matrix  and 
1x10  11  m2  which  is  the  fracture  zone  permeability  for  the  three  fractures. 
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D.5  Mole  Fraction  of  TCE 


To  include  TCE  in  the  model,  TMVOC  requires  it  to  be  in  the  form  of  a  mole  fraction.  The  initial 
concentration  of  TCE  estimated  from  the  pre-heating  rock  samples  was  used  to  calculate  the  values  used 
in  the  model.  It  was  assumed  that  only  the  sorbed  and  aqueous  phase  is  present  and  the  following 
equation  was  used  to  calculate  the  aqueous  TCE  concentration  (Kueper  and  Davies,  2009): 


Ct  =  ^~  (KdPb  +  0w  +  H'0a)  (D‘9) 

Pb 

where  Cris  the  soil  concentration,  the  dry  soil  bulk  density  of  2520  kg/m3,  Kd  is  the  soil-water  partition 
coefficient  and  is  equal  to  koc  ■  foc.  The  value  used  for  the  organic  carbon  partition  coefficient,  koc,  is  126 
mg/L  which  is  the  standard  value  for  TCE,  and  foc,  the  fraction  of  organic  carbon  measured  by  the 
laboratory  was  found  to  be  0.00785.  The  water-filled  porosity,  0W,  is  equal  to  the  porosity  as  the  soil  is 
assumed  to  be  saturated  and  is  0.033.  Since  there  is  no  air  present  in  the  system,  the  air-filled  porosity, 
0a,  is  set  to  zero.  Rearranging  Equation  D.9  to  solve  for  C„  the  effective  solubility  in  mg/L  is  as  follows: 


CrPi) 

KdPb  +  0w 


(D.10) 
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The  effective  solubility  can  be  converted  into  units  of  mol/L  by  dividing  C,  by  the  molecular  weight  of 
trichloroethylene  (131.39  g/mol).  The  mole  fraction  is  calculated  using  Equation  D.ll  below  and  a 
molecular  weight  of  water  (CM)  of  55.46  mol/L. 


[TCE]  = 


Cj 

Cw  +  Ci 


(D.ll) 
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Appendix  E 


Fixed  versus  No-Flow  Bottom  Boundary  Condition 

The  bottom  boundary  in  numerical  model  simulations  was  extended  10.5  m  below  the  bottom  of  the 
treatment  zone  in  order  to  reduce  the  influence  of  the  bottom  boundary  condition  on  the  model  results. 
There  were  two  different  bottom  boundary  conditions  examined  when  setting  up  the  baseline  model 
simulations,  a  fixed  boundary  and  a  no-flow  boundary. 

A  fixed  boundary  is  set  with  a  constant  pressure  and  temperature  during  the  model  simulations  which 
forces  those  boundary  conditions  to  always  be  met  at  that  boundary.  A  no-flow  boundary  is  essentially  a 
confining  layer  that  does  not  allow  any  fluid  or  temperature  cross  the  boundaiy.  The  baseline  model  (with 
a  matrix  permeability  of  1  x  1 0  1  m2)  was  run  with  both  bottom  boundary  conditions  and  compared  to 
examine  the  effects  the  different  boundaries  created  on  the  TCE  removal  from  the  rock  matrix. 

The  two  model  simulations  were  compared  at  the  observation  location  for  Figures  E-l  to  E-3  to  examine 
the  difference  in  temperature,  pressure,  gas  saturation  and  TCE  extraction  at  four  different  time  steps: 

•  30  days  which  is  in  the  heating  phase; 

•  40  days  which  is  when  the  maximum  average  temperature  is  reached;  and, 

•  50  and  80  days  during  the  period  of  constant  heating. 

Figure  E-l  demonstrates  that  there  is  no  significant  difference  in  the  temperature  profiles  at  the 
observation  location  between  the  two  boundary  conditions. 
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30  Days 


40  Days 


Fixed  Bottom  Boundary  Condition 
No  Flow  Bottom  Boundary  Condition 
Fracture 

Bottom  of  Heated  Zone 


Figure  E-l  -  Model-predicted  temperature  profile  in  the  rock  matrix  at  the  observation  location 

over  time. 


A  similarity  in  pressure  profiles  was  also  seen  at  the  observation  location  between  the  two  boundary 
conditions  (Figure  E-2). 
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Fixed  Bottom  Boundary  Condition 
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Figure  E-2  -  Model-predicted  pressure  profile  in  the  rock  matrix  at  the  observation  location  over 

time. 


The  primary  difference  was  displayed  through  the  gas  saturation  profdes.  Figure  E-3  displays  the 
difference  between  the  gas  saturation  profiles  in  the  fixed  bottom  boundary  and  no-flow  bottom  boundary 
simulations  at  the  observation  location.  The  figure  displays  that  the  boundaries  have  a  significant  effect 
on  the  gas  saturation  in  the  rock  matrix  during  heating. 
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30  Days  40  Days 


-  Fixed  Bottom  Boundary  Condition 

- No-Flow  Bottom  Boundary  Condition 

Fracture 


Figure  E-3  -  Model-predicted  gas  saturation  profile  in  the  rock  matrix  at  the  observation  location 

over  time. 


The  profiles  show  that  at  30  days  there  is  little  difference  between  the  gas  saturation  profiles  with  the 
different  boundary  conditions.  When  the  maximum  temperature  is  reached  at  40  days,  the  saturation 
profile  becomes  quite  different  between  the  two  boundaries.  The  fixed  bottom  boundary  condition  shows 


a  less  gas  saturation  in  the  rock  matrix  than  the  no-flow  bottom  boundary  condition  for  40,  50  and  80 
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days.  This  is  because  the  fixed  boundary  condition  allows  the  model  to  “fill”  back  up  with  water  therefore 
reducing  the  gas  phase  (this  is  due  to  the  constant  pressure  being  fixed  at  the  bottom  boundary). 


When  the  TCE  in  the  rock  matrix  is  converted  from  the  liquid  phase  to  the  gas  phase  it  is  easier  to  extract. 
Figure  E-4  displays  the  TCE  profile  in  the  rock  matrix  over  time.  The  TCE  concentration  for  the  no-flow 
bottom  boundaiy  is  always  less  than  that  observed  in  the  fixed  bottom  boundaiy  simulation.  The  no-flow 
bottom  boundary  resulted  in  67%  TCE  removal  from  the  rock  matrix  compared  to  the  fixed  bottom 
boundaiy  which  had  44%  TCE  removal. 
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30  Days  40  Days 


-  Fixed  Bottom  Boundary  Condition 

- No-Flow  Bottom  Boundary  Condition 

Fracture 


Figure  E-4  -  Model-predicted  TCE  profile  in  the  rock  matrix  at  the  observation  location  over  time. 


Although  both  boundary  conditions  are  located  10.5  m  away  from  the  heating  zone  there  is  a  significant 
effect  on  the  amount  of  gas  that  is  created  in  the  rock  matrix  during  heating.  This  difference  allows  for 
more  TCE  to  be  removed  from  the  rock  matrix  in  the  no-flow  bottom  boundary  than  the  fixed  bottom 
boundary.  The  average  percentage  of  TCE  removed  from  the  rock  matrix  during  the  NAWC  pilot  test 


after  100  days  of  heating  was  63.5%.  The  rock  matrix  at  the  NAWC  site  becomes  less  fractured  with 

98 


depth  therefore  it  could  be  hypothesized  that  the  lower  permeability  with  depth  acts  to  restrict  the  upward 
flow  of  water  back  into  the  heating  zone.  The  no-flow  bottom  boundary  condition  acts  in  a  similar  way  as 
a  lower  permeability  zone  which  resulted  in  67%  removal  from  the  rock  matrix.  These  results  concluded 
that  a  no-flow  bottom  boundary  is  more  appropriate  to  use  than  the  fixed  bottom  boundary  condition  and 
was  therefore  applied  for  all  simulations. 


99 


